Associations Between Neighbourhood Socioeconomic Status and Road Traffic Injury Risk in Urban Areas of British Columbia - R Code
Michael Branion-Calles
2024-10-07
Load libraries and functions
# install.packages("INLA",repos=c(getOption("repos"),INLA="https://inla.r-inla-download.org/R/stable"), dep=TRUE)
start_time <- Sys.time()
library(INLA) # For fitting integrated nested Laplace approximation (INLA) models
library(tidyverse) # For data manipulation and visualization
library(sp) # For spatial data handling
library(spdep) # For spatial dependence analysis
library(flextable) # For creating flexible tables
library(RColorBrewer) # For color palettes
library(rcartocolor) # For additional color palettes
library(cowplot) # For combining plots
library(janitor) # For cleaning data
library(broom)
# Functions ---------------------------------------------------------------
# Standardize a numeric vector
my_std <- function(x) { (x - mean(x,na.rm=TRUE)) / sd(x,na.rm = TRUE)}
# Fit a BYM2 model using INLA
fit_bym2 <- function(table,
outcome_var,
interest_var,
distribution = "poisson",
formula_terms = NULL,
prec_param = c(1, 0.001),
phi_param = c(0.5, 0.5)) {
# Define priors for precision and phi
prior <- list(
prec = list(
prior = "pc.prec",
param = prec_param),
phi = list(
prior = "pc",
param = phi_param)
)
# Construct the INLA formula based on input parameters
if (is.null(formula_terms)) {
inla_formula <- as.formula(paste(outcome_var, "~", interest_var, "+ f(ui, model = 'bym2', graph = g, hyper = prior, adjust.for.con.comp = TRUE,
constr = TRUE, scale.model = TRUE)"))
} else {
inla_formula <- as.formula(paste(outcome_var, "~", interest_var, "+", formula_terms, "+ f(ui, model = 'bym2', graph = g, hyper = prior, adjust.for.con.comp = TRUE,
constr = TRUE, scale.model = TRUE)"))
}
#adjust.for.con.comp = TRUE, scale.model = TRUE: assigns one intercept to each region in addition to using a sum-to-zero constraint for each connected region. By adding an intercept for each region, we infer that the baseline prevalence is different in the disconnected regions.
#The spatial random effects can be interpreted as the area-specific deviation from the region-specific risk in this case. The intercepts for the disconnected regions need to be explicitly specified in INLA:
# Fit the INLA model with BYM2 prior
fit <- inla(inla_formula,
data = table,
family = distribution,
control.compute = list(dic = TRUE, waic = TRUE,cpo=TRUE,config=TRUE),
control.predictor = list(compute = TRUE),
verbose = FALSE)
# Return the INLA model object
return(fit)
}
# Tidy up INLA model fixed effects results
tidy_inla <- function(fixed_effects, exponentiate = FALSE, sigfig =NULL) {
# Create a tidy data frame
tidy_results <- fixed_effects %>%
dplyr::rename(estimate = 'mean',
std.error = 'sd',
lower_ci = '0.025quant',
upper_ci = '0.975quant')
if (exponentiate) {
tidy_results <- tidy_results %>%
mutate(estimate = exp(estimate),
lower_ci = exp(lower_ci),
upper_ci = exp(upper_ci))
}
if (!is.null(sigfig)) {
tidy_results <- tidy_results %>%
mutate(estimate = signif(estimate, sigfig),
std.error = signif(std.error, sigfig),
lower_ci = signif(lower_ci, sigfig),
upper_ci = signif(upper_ci, sigfig))
}
return(tidy_results)
}
# prs <- function(x,n.dig = 2){
#
# formatC(round(x,n.dig), format='f', big.mark = ",",digits = n.dig)
#
# }
cma_models <- function(cma_data, cma_name, outcome, deprivation_indicator,distribution = "poisson") {
###### Spatial Models #####
### 1: Unadjusted
print("Fit unadjusted")
m1 <- fit_bym2(cma_data, outcome, deprivation_indicator,distribution)
m1_fe <- tidy_inla(m1$summary.fixed, exponentiate = TRUE, sigfig =3) %>%
mutate(
`IRR (95% CI)` = paste0(estimate, " (", lower_ci, ", ", upper_ci, ")"),
variable = row.names(.)
) %>%
select(variable, `IRR (95% CI)`)
print("Fit adjusted by road length")
### 2: Adjusted by Road Length
m2 <- fit_bym2(cma_data, outcome, deprivation_indicator, formula_terms = "ln_roads_km_c",distribution)
m2_fe <- tidy_inla(m2$summary.fixed, exponentiate = TRUE, sigfig =3) %>%
mutate(
`IRR (95% CI)` = paste0(estimate, " (", lower_ci, ", ", upper_ci, ")"),
variable = row.names(.)
) %>%
select(variable, `IRR (95% CI)`)
### 3: Adjusted by Road Length and other variables
m3 <- fit_bym2(cma_data, outcome, deprivation_indicator,
formula_terms = "ln_roads_km_c + roads_prop_highway_arterial_z + ale_index_z + canbics_index_z + population_100_c",distribution)
m3_fe <- tidy_inla(m3$summary.fixed, exponentiate = TRUE, sigfig =3) %>%
mutate(
`IRR (95% CI)` = paste0(estimate, " (", lower_ci, ", ", upper_ci, ")"),
variable = row.names(.)
) %>%
select(variable, `IRR (95% CI)`)
l <- list(m1, m2, m3)
names(l) <- c("Unadjusted","Adjusted1","Adjusted2")
print("Fit adjusted by road length and covariates")
# Compile results into a tibble
tibble(
cma = cma_name,
model_type = c("unadjusted", "adjusted_minimal", "adjusted_full"),
models = l,
fixed_effects = list(m1_fe, m2_fe, m3_fe)
)
}
# Clean and format fixed effects results from models
clean_fixed_effects <- function(cma_model_object){
cma_model_object %>%
bind_rows(,.id="model_type") %>%
mutate(model_type = case_when(model_type == "1"~"Unadjusted IRR (95% CI)",
model_type == "2"~"Minimally Adjusted IRR (95% CI)",
model_type == "3"~"Fully Adjusted IRR (95% CI)",
)) %>%
pivot_wider(id_cols = variable,names_from = model_type,values_from = `IRR (95% CI)` ) %>%
mutate(variable = case_when(variable == "ln_roads_km_c" ~"Log(Total KMs of Road)",
variable == "roads_prop_highway_arterial_z" ~ "% of roads classifed as Arterial/Highway (SD)",
variable == "ale_index_z" ~ "Canadian Active Living Enrivonment Index (SD)",
variable == "canbics_index_z" ~ "Canadian Bikeway Safety and Comfort Index (SD)",
variable == "population_100_c" ~ "Population (100)",
TRUE ~ variable
))
}
assign_nearest_neighbors <- function(nb_object, sp_object,make_symmetrical=TRUE) {
# Calculate centroids for isolated polygons
centroids <- coordinates(sp_object) # If using sp
centroids_nn <- knearneigh(centroids, k = 1)
# Get the list of neighbors
zero_nb_logical <- sapply(nb_object, function(x) x[1] ==0)
zero_nb_indices <- which(zero_nb_logical)
zero_nn <- centroids_nn$nn[zero_nb_indices, ]
for(i in 1:length(zero_nn)){
nb_object[zero_nb_indices[i]] <- zero_nn[i]
}
if(make_symmetrical==TRUE){
nb_object <- make.sym.nb(nb_object)
}
return(nb_object)
}
create_descriptive_table <- function(sf_obj, selected_vars) {
sf_obj %>%
# Drop geometry and select only the specified variables
st_drop_geometry() %>%
select(all_of(selected_vars)) %>%
# Gather descriptive statistics for each variable
summarise(across(everything(), list(
sum = ~sum(.x, na.rm = TRUE),
mean = ~mean(.x, na.rm = TRUE),
sd = ~sd(.x, na.rm = TRUE),
min = ~min(.x, na.rm = TRUE),
max = ~max(.x, na.rm = TRUE),
n = ~sum(!is.na(.x)|is.na(.x)), # Count non-NA values
missing = ~sum(is.na(.x)),
n_zero = ~sum(.x==0,na.rm=TRUE)# Count 0 values
), .names = "{.col}-{.fn}")) %>%
# Pivot longer to have each variable as its own row
pivot_longer(everything(),
names_to = c("variable", ".value"),
names_sep = "-") %>%
# Add a column to identify the sf object
# Reorder columns for better readability
select(variable, n, sum, mean, sd, min, max, , missing,n_zero)
}Load Data and Clean
# Load and clean data -----------------------------------------------------
cma_name <- cancensus::get_census("CA21",
level = "CMA",
regions = list(PR=59)) %>%
rename(CMA_UID = GeoUID,
`CMA Name` = `Region Name`,
cma_population = Population) %>%
arrange(desc(cma_population)) %>%
select(CMA_UID, Type, `CMA Name`,cma_population) %>%
mutate(`CMA Name` = str_remove_all(`CMA Name`," \\(B\\)|\\(K\\)|\\(D\\)")) %>%
clean_names() %>%
mutate(cma_name = stringr::str_trim(cma_name,"right"))## Downloading: 920 B Downloading: 920 B Downloading: 920 B Downloading: 920 B Downloading: 920 B Downloading: 920 B# Read spatial data for dissemination areas in British Columbia
bc_da <- st_read(
"C:/Users/micha/Documents/GitHub/Area-Level-Deprivation-Traffic-Injury/Processed Data/da_v3_2021.gpkg"
) %>%
# Join the census data with the spatial data
left_join(cma_name, by = join_by(cma_uid)) %>%
filter(!is.na(cma_name)) %>%
mutate(
population_100 = population/100,
population_100_c = scale(population_100, scale = FALSE), # Convert road lengths from meters to kilometers
total_roads_km = total_roads_m / 1000,
n_casualty_claims_rate = n_casualty_claims/total_roads_km*10,
n_pedestrian_casualty_claims_rate = n_pedestrian_casualty_claims/total_roads_km*10,
n_cyclist_casualty_claims_rate = n_cyclist_casualty_claims/total_roads_km*10,
# Calculate the proportion of roads that are highways or arterials
roads_prop_highway_arterial = (highway_m + arterial_m) / total_roads_m,
# Standardize the proportion of highways and arterials
roads_prop_highway_arterial_z = my_std(roads_prop_highway_arterial),
# Log-transform the total length of roads in kilometers
ln_roads_km = ifelse(total_roads_km ==0, NA,log(total_roads_km)),
ln_roads_km_c = scale(ln_roads_km, scale = FALSE),
# Calculate the casualty claims rate by total road kilometers
casualty_claims_rate = n_casualty_claims / total_roads_km * 10,
# Calculate the cyclist casualty claims rate by total road kilometers
cyclist_casualty_claims_rate = n_cyclist_casualty_claims / total_roads_km * 10,
# Calculate the pedestrian casualty claims rate by total road kilometers
pedestrian_casualty_claims_rate = n_pedestrian_casualty_claims / total_roads_km * 10,
# Convert journey-to-work total from meters to kilometers
jtw_total_100 = jtw_total / 100,
# Standardize the journey-to-work total
jtw_total_z = my_std(jtw_total),
# Calculate the proportion of journey-to-work trips that are motor vehicle-related
jtw_prop_mv = ((jtw_driver + jtw_passenger) / jtw_total),
# Calculate the proportion of non-motor vehicle journey-to-work trips
jtw_prop_non_mv = (1 - jtw_prop_mv),
# Standardize the proportion of non-motor vehicle trips
jtw_prop_non_mv_z = my_std(jtw_prop_non_mv),
# Adjust and scale the remote index for visualization
remote_index_2 = (1 - remote_index) * 100,
# Log-transform the adjusted remote index
ln_remote_index = log(remote_index_2),
# Standardize the VANDIX index
vandix_z = my_std(vandix),
# Convert VANDIX quintile to descriptive categories
vandix_quintile = vandix_quintile %>% as.character(),
vandix_quintile = case_when(
vandix_quintile == "1" ~ "1 - Least Deprived",
vandix_quintile == "5" ~ "5 - Most Deprived",
TRUE ~ vandix_quintile
),
# Convert educational attainment quintile to descriptive categories
no_highschool_quintile = no_highschool_quintile %>% as.character(),
no_highschool_quintile = case_when(
no_highschool_quintile == "1" ~ "1 - Least Deprived",
no_highschool_quintile == "5" ~ "5 - Most Deprived",
TRUE ~ no_highschool_quintile
),
# Convert unemployment rate quintile to descriptive categories
unemployment_rate_quintile = unemployment_rate_quintile %>% as.character(),
unemployment_rate_quintile = case_when(
unemployment_rate_quintile == "1" ~ "1 - Least Deprived",
unemployment_rate_quintile == "5" ~ "5 - Most Deprived",
TRUE ~ unemployment_rate_quintile
),
# Convert university degree quintile to descriptive categories
university_degree_quintile = university_degree_quintile %>% as.character(),
university_degree_quintile = case_when(
university_degree_quintile == "1" ~ "1 - Least Deprived",
university_degree_quintile == "5" ~ "5 - Most Deprived",
TRUE ~ university_degree_quintile
),
# Convert lone-parent family quintile to descriptive categories
lone_parent_fam_quintile = lone_parent_fam_quintile %>% as.character(),
lone_parent_fam_quintile = case_when(
lone_parent_fam_quintile == "1" ~ "1 - Least Deprived",
lone_parent_fam_quintile == "5" ~ "5 - Most Deprived",
TRUE ~ lone_parent_fam_quintile
),
# Convert average household income quintile to descriptive categories
hh_avg_income_quintile = hh_avg_income_quintile %>% as.character(),
hh_avg_income_quintile = case_when(
hh_avg_income_quintile == "1" ~ "1 - Least Deprived",
hh_avg_income_quintile == "5" ~ "5 - Most Deprived",
TRUE ~ hh_avg_income_quintile
),
# Convert home ownership quintile to descriptive categories
home_owner_quintile = home_owner_quintile %>% as.character(),
home_owner_quintile = case_when(
home_owner_quintile == "1" ~ "1 - Least Deprived",
home_owner_quintile == "5" ~ "5 - Most Deprived",
TRUE ~ home_owner_quintile
),
# Convert participation rate quintile to descriptive categories
participation_rate_quintile = participation_rate_quintile %>% as.character(),
participation_rate_quintile = case_when(
participation_rate_quintile == "1" ~ "1 - Least Deprived",
participation_rate_quintile == "5" ~ "5 - Most Deprived",
TRUE ~ participation_rate_quintile
),
# Convert CANBICS class to descriptive categories
canbics_class_c = case_when(
canbics_class == 1 | canbics_class == "2" ~ "1 - Low",
canbics_class == 5 | canbics_class == "4" ~ "3 - High",
canbics_class == 3 ~ "2 - Medium"
),
# Standardize the CANBICS index
canbics_index_z = my_std(canbics_index),
# Convert ALE class to descriptive categories
ale_class_c = case_when(
ale_class == 1 | ale_class == "2" ~ "1 - Low",
ale_class == 5 | ale_class == "4" ~ "3 - High",
ale_class == 3 ~ "2 - Medium"
),
# Standardize the ALE index
ale_index_z = my_std(ale_index),
# Calculate VanDIX
z_no_highschool_prevalance_w = scale(no_highschool_prevalance),
z_university_degree_prevalance_w = scale(university_degree_prevalance * -1),
z_unemployment_rate_w = scale(unemployment_rate),
z_lone_parent_fam_prevalence_w = scale(lone_parent_fam_prevalence),
z_hh_avg_income_w = scale(hh_avg_income * -1),
z_home_owner_prevalence_w = scale(home_owner_prevalence * -1),
z_participation_rate_w = scale(participation_rate * -1),
# Compute VanDIX (Area Level Deprivation Index) score
vandix = as.numeric(
z_hh_avg_income_w * 0.089 +
z_home_owner_prevalence_w * 0.089 +
z_lone_parent_fam_prevalence_w * 0.143 +
z_no_highschool_prevalance_w * 0.25 +
z_university_degree_prevalance_w * 0.179 +
z_unemployment_rate_w * 0.214 +
z_participation_rate_w * 0.036
),
vandix_z_c = case_when(vandix_z< -5 ~ "<-5",
vandix_z>= -5 & vandix_z< -4 ~"-5 - -4",
vandix_z>= -4 & vandix_z < -3 ~"-4 - -3",
vandix_z>= -3 & vandix_z < -2 ~"-3 - -2",
vandix_z>= -2 & vandix_z < -1 ~"-2 - -1",
vandix_z>= -1 & vandix_z < 0 ~"-1 - 0",
vandix_z>= 0 & vandix_z < 1 ~"0 - 1",
vandix_z>= 1 & vandix_z<2 ~"1 - 2",
vandix_z>= 2 & vandix_z<3 ~"2 - 3",
vandix_z>= 3 & vandix_z<4 ~"3 - 4",
vandix_z>= 4 & vandix_z<5 ~"4 - 5",
vandix_z >= 5 ~ ">5",
),
vandix_z_c = factor(vandix_z_c,levels = c("<-5","-5 - -4","-4 - -3" , "-3 - -2" , "-2 - -1" , "-1 - 0", "0 - 1" , "1 - 2", "2 - 3", "3 - 4","4 - 5",">5"))
) %>%
mutate(region_name = case_when(
cma_name %in% c("Prince Rupert", "Terrace") ~ "Northwest",
cma_name %in% c("Fort St. John", "Dawson Creek", "Prince George", "Quesnel", "Williams Lake") ~ "North Central",
cma_name %in% c("Kamloops", "Salmon Arm") ~ "Kamloops-Salmon Arm",
cma_name %in% c("Vernon", "Kelowna", "Penticton") ~ "Okanagan",
cma_name %in% c("Cranbrook", "Nelson", "Trail") ~ "Southeast",
cma_name %in% c("Vancouver", "Squamish") ~ "Vancouver-Squamish",
cma_name %in% c("Abbotsford - Mission", "Chilliwack") ~ "Fraser Valley",
cma_name %in% c("Campbell River", "Courtenay", "Port Alberni", "Parksville", "Nanaimo", "Ladysmith", "Duncan","Powell River") ~ "Central Island-Powell River",
cma_name == "Victoria" ~ "Victoria",
TRUE ~ NA_character_ # Handle cases where cma_name doesn't match any of the specified values
))## Reading layer `da_v3_2021' from data source
## `C:\Users\micha\Documents\GitHub\Area-Level-Deprivation-Traffic-Injury\Processed Data\da_v3_2021.gpkg'
## using driver `GPKG'
## Simple feature collection with 7848 features and 130 fields
## Geometry type: MULTIPOLYGON
## Dimension: XY
## Bounding box: xmin: 273876.1 ymin: 367780.7 xmax: 1870608 ymax: 1735671
## Projected CRS: NAD83 / BC Albers#split data into list of CMAs
bc_cma <- bc_da %>%
split(.$cma_name)
bc_region <- bc_da %>%
split(.$region_name)
cma_pop <- bc_da %>%
st_drop_geometry() %>%
group_by(cma_name) %>%
summarise(n_da = n(),
population = sum(population,na.rm=TRUE)) %>%
arrange(desc(population))
region_pop <- bc_da %>%
st_drop_geometry() %>%
group_by(region_name) %>%
summarise(n_da = n(),
population = sum(population,na.rm=TRUE)) %>%
arrange(desc(population))
region_descriptives <- map(
bc_region,
create_descriptive_table,
c(
"n_claims",
"n_casualty_claims",
"n_cyclist_claims",
"n_cyclist_casualty_claims",
"n_pedestrian_claims",
"n_pedestrian_casualty_claims",
"population",
"total_roads_km",
"roads_prop_highway_arterial",
"no_highschool_prevalance",
"unemployment_rate",
"hh_avg_income",
"participation_rate",
"university_degree_prevalance",
"lone_parent_fam_prevalence",
"home_owner_prevalence",
"vandix",
"roads_prop_highway_arterial",
"ale_index",
"canbics_index"
)
) %>%
bind_rows(.id = "region") %>%
arrange(desc(n)) %>%
mutate(p_zero = n_zero/n*100)
cma_descriptives <- map(
bc_cma,
create_descriptive_table,
c(
"n_claims",
"n_casualty_claims",
"n_cyclist_claims",
"n_cyclist_casualty_claims",
"n_pedestrian_claims",
"n_pedestrian_casualty_claims",
"population",
"total_roads_km",
"roads_prop_highway_arterial",
"no_highschool_prevalance",
"unemployment_rate",
"hh_avg_income",
"participation_rate",
"university_degree_prevalance",
"lone_parent_fam_prevalence",
"home_owner_prevalence",
"vandix",
"roads_prop_highway_arterial",
"ale_index",
"canbics_index"
)
) %>%
bind_rows(.id = "cma") %>%
arrange(desc(n)) %>%
mutate(p_zero = n_zero/n*100)Spatial Models
Vancouver - Squamish
Descriptive Statistics
region | variable | n | sum | mean | sd | min | max | missing | n_zero | p_zero |
|---|---|---|---|---|---|---|---|---|---|---|
Vancouver-Squamish | n_claims | 3,622 | 693,530.0 | 191.5 | 363.2 | 0.0 | 6,315.0 | 0 | 12 | 0.3 |
n_casualty_claims | 149,692.0 | 41.3 | 84.8 | 0.0 | 1,896.0 | 0 | 133 | 3.7 | ||
n_cyclist_claims | 7,656.0 | 2.1 | 4.7 | 0.0 | 115.0 | 0 | 1,457 | 40.2 | ||
n_cyclist_casualty_claims | 4,968.0 | 1.4 | 3.2 | 0.0 | 72.0 | 0 | 1,859 | 51.3 | ||
n_pedestrian_claims | 9,539.0 | 2.6 | 5.1 | 0.0 | 80.0 | 0 | 1,290 | 35.6 | ||
n_pedestrian_casualty_claims | 7,502.0 | 2.1 | 4.0 | 0.0 | 68.0 | 0 | 1,458 | 40.3 | ||
population | 2,667,057.0 | 736.6 | 532.3 | 0.0 | 8,800.0 | 1 | 10 | 0.3 | ||
total_roads_km | 15,056.8 | 4.2 | 5.7 | 0.0 | 176.7 | 0 | 2 | 0.1 | ||
roads_prop_highway_arterial | 577.0 | 0.2 | 0.2 | 0.0 | 1.0 | 2 | 1,426 | 39.4 | ||
no_highschool_prevalance | 476.5 | 0.1 | 0.1 | 0.0 | 0.8 | 300 | 16 | 0.4 | ||
unemployment_rate | 19,517.2 | 5.9 | 3.6 | 0.0 | 40.0 | 300 | 339 | 9.4 | ||
hh_avg_income | 279,328,113.0 | 85,421.4 | 34,171.3 | 0.0 | 560,267.0 | 352 | 3 | 0.1 | ||
participation_rate | 215,750.1 | 64.9 | 10.1 | 0.0 | 94.6 | 300 | 1 | 0.0 | ||
university_degree_prevalance | 1,890.5 | 0.6 | 0.1 | 0.1 | 1.0 | 300 | 0 | 0.0 | ||
lone_parent_fam_prevalence | 513.7 | 0.2 | 0.1 | 0.0 | 0.7 | 299 | 15 | 0.4 | ||
home_owner_prevalence | 2,229.1 | 0.7 | 0.2 | 0.0 | 1.3 | 300 | 36 | 1.0 | ||
vandix | -473.2 | -0.1 | 0.6 | -2.1 | 3.9 | 352 | 0 | 0.0 | ||
ale_index | 5,407.4 | 1.6 | 3.6 | -2.1 | 25.1 | 310 | 13 | 0.4 | ||
canbics_index | 15,159.2 | 4.6 | 4.1 | 0.0 | 24.3 | 305 | 72 | 2.0 |
Maps
claims <- ggplot() +
geom_sf(data = vancouver_sf, aes(fill = n_casualty_claims,colour=n_casualty_claims)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
theme_void() +
ggtitle(
"Vancouver"
)
vandix <- ggplot() +
geom_sf(data = vancouver_sf, aes(fill = vandix_z_c,colour=vandix_z_c)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
scale_colour_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
theme_void()
total_roads <- ggplot() +
geom_sf(data = vancouver_sf, aes(fill = total_roads_km,colour=total_roads_km)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
theme_void()
cowplot::plot_grid(claims,vandix,total_roads,ncol=1)
Setting up Spatial Weights Matrix
#### Define Spatial Neighrbourhoods
vancouver_sp <- as(vancouver_sf, "Spatial")
vancouver_sp$ui <- 1:nrow(vancouver_sp@data)
coords <- coordinates(vancouver_sp)
vancouver_nb <- poly2nb(vancouver_sp, queen = TRUE)## Warning in poly2nb(vancouver_sp, queen = TRUE): some observations have no neighbours;
## if this seems unexpected, try increasing the snap argument.## Warning in poly2nb(vancouver_sp, queen = TRUE): neighbour object has 6 sub-graphs;
## if this sub-graph count seems unexpected, try increasing the snap argument.vancouver_nb## Neighbour list object:
## Number of regions: 3622
## Number of nonzero links: 22972
## Percentage nonzero weights: 0.1751064
## Average number of links: 6.342352
## 1 region with no links:
## 2938
## 6 disjoint connected subgraphs#assign nearest neighbour for no links
vancouver_nb <- assign_nearest_neighbors(vancouver_nb,vancouver_sp)
plot(vancouver_sp, border = grey(0.5))
plot(vancouver_nb,
coords = coords,
add = TRUE, pch = 16, lwd = 2)
Spatial Dependency in each crash type
listw <- nb2listw(vancouver_nb,zero.policy = TRUE)
all_cc_mi <- moran.test(vancouver_sf$n_casualty_claims, listw,zero.policy = TRUE)
all_cc_mi##
## Moran I test under randomisation
##
## data: vancouver_sf$n_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 24.441, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 2.318516e-01 -2.761668e-04 9.020552e-05cyc_cc_mi <- moran.test(vancouver_sf$n_cyclist_casualty_claims, listw,zero.policy = TRUE)
cyc_cc_mi##
## Moran I test under randomisation
##
## data: vancouver_sf$n_cyclist_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 38.928, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 3.683422e-01 -2.761668e-04 8.966834e-05pd_cc_mi <- moran.test(vancouver_sf$n_pedestrian_casualty_claims, listw,zero.policy = TRUE)
pd_cc_mi##
## Moran I test under randomisation
##
## data: vancouver_sf$n_pedestrian_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 32.695, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 3.122914e-01 -2.761668e-04 9.139466e-05BYM2 Models
nb2INLA("vancouver.adj", vancouver_nb)
g <- inla.read.graph(filename = "vancouver.adj")
# Model Set 1: Total Casualty Claim Crashes
van_models_1 <- cma_models(vancouver_sp@data, "Vancouver-Squamish", "n_casualty_claims", "vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"van_all <- clean_fixed_effects(van_models_1$fixed_effects)
map(van_models_1$models, ~ summary(.x))## $Unadjusted
## Time used:
## Pre = 18.7, Running = 5.27, Post = 0.528, Total = 24.5
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 2.802 0.016 2.771 2.802 2.834 2.802 0
## vandix_z 0.211 0.033 0.147 0.211 0.276 0.211 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.311 0.014 0.284 0.310 0.339 0.310
## Phi for ui 0.794 0.021 0.751 0.795 0.833 0.797
##
## Deviance Information Criterion (DIC) ...............: 23817.67
## Deviance Information Criterion (DIC, saturated) ....: 7368.64
## Effective number of parameters .....................: 3420.00
##
## Watanabe-Akaike information criterion (WAIC) ...: 23712.55
## Effective number of parameters .................: 2288.45
##
## Marginal log-Likelihood: -15034.43
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 18.6, Running = 4.79, Post = 0.636, Total = 24
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 3.254 0.015 3.223 3.254 3.284 3.254 0
## vandix_z 0.327 0.026 0.276 0.327 0.378 0.327 0
## ln_roads_km_c 1.294 0.029 1.238 1.294 1.351 1.294 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.508 0.020 0.470 0.508 0.548 0.507
## Phi for ui 0.811 0.016 0.778 0.811 0.841 0.812
##
## Deviance Information Criterion (DIC) ...............: 23705.93
## Deviance Information Criterion (DIC, saturated) ....: 7256.88
## Effective number of parameters .....................: 3264.92
##
## Watanabe-Akaike information criterion (WAIC) ...: 23514.43
## Effective number of parameters .................: 2154.47
##
## Marginal log-Likelihood: -14241.75
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 18.2, Running = 5.04, Post = 0.736, Total = 24
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode
## (Intercept) 2.989 0.017 2.957 2.989 3.022 2.989
## vandix_z 0.222 0.022 0.180 0.222 0.265 0.222
## ln_roads_km_c 0.990 0.027 0.937 0.990 1.043 0.990
## roads_prop_highway_arterial_z 0.593 0.016 0.562 0.593 0.624 0.593
## ale_index_z 0.109 0.037 0.036 0.109 0.182 0.109
## canbics_index_z 0.209 0.037 0.137 0.209 0.282 0.209
## population_100_c 0.017 0.003 0.012 0.017 0.022 0.017
## kld
## (Intercept) 0
## vandix_z 0
## ln_roads_km_c 0
## roads_prop_highway_arterial_z 0
## ale_index_z 0
## canbics_index_z 0
## population_100_c 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.838 0.038 0.765 0.837 0.915 0.835
## Phi for ui 0.777 0.022 0.733 0.778 0.818 0.780
##
## Deviance Information Criterion (DIC) ...............: 23560.95
## Deviance Information Criterion (DIC, saturated) ....: 7111.89
## Effective number of parameters .....................: 3108.88
##
## Watanabe-Akaike information criterion (WAIC) ...: 23296.41
## Effective number of parameters .................: 2021.13
##
## Marginal log-Likelihood: -13603.57
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 2: Total Casualty Cyclist Claim Crashes
van_models_2 <- cma_models(vancouver_sp@data,"Vancouver-Squamish","n_cyclist_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"van_cyclist <- clean_fixed_effects(van_models_2$fixed_effects)
map(van_models_2$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 18.3, Running = 3.99, Post = 0.406, Total = 22.6
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -0.523 0.032 -0.588 -0.523 -0.461 -0.523 0
## vandix_z 0.148 0.037 0.075 0.148 0.221 0.148 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.567 0.042 0.489 0.565 0.653 0.562
## Phi for ui 0.633 0.048 0.534 0.635 0.724 0.638
##
## Deviance Information Criterion (DIC) ...............: 9767.72
## Deviance Information Criterion (DIC, saturated) ....: 5213.94
## Effective number of parameters .....................: 1755.44
##
## Watanabe-Akaike information criterion (WAIC) ...: 9992.02
## Effective number of parameters .................: 1391.74
##
## Marginal log-Likelihood: -4164.36
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 18.2, Running = 4.69, Post = 0.463, Total = 23.4
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -0.134 0.026 -0.186 -0.134 -0.083 -0.134 0
## vandix_z 0.217 0.032 0.154 0.217 0.280 0.217 0
## ln_roads_km_c 1.072 0.035 1.003 1.072 1.141 1.072 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.622 0.047 0.532 0.621 0.718 0.621
## Phi for ui 0.934 0.026 0.874 0.938 0.974 0.945
##
## Deviance Information Criterion (DIC) ...............: 9035.13
## Deviance Information Criterion (DIC, saturated) ....: 4481.34
## Effective number of parameters .....................: 1133.23
##
## Watanabe-Akaike information criterion (WAIC) ...: 9017.22
## Effective number of parameters .................: 857.84
##
## Marginal log-Likelihood: -3785.23
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 18.4, Running = 5.09, Post = 0.709, Total = 24.2
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -0.291 0.031 -0.352 -0.291 -0.231
## vandix_z 0.185 0.032 0.123 0.185 0.248
## ln_roads_km_c 0.844 0.039 0.767 0.844 0.921
## roads_prop_highway_arterial_z 0.301 0.024 0.254 0.301 0.348
## ale_index_z 0.041 0.045 -0.048 0.041 0.130
## canbics_index_z 0.092 0.048 -0.001 0.092 0.186
## population_100_c 0.022 0.003 0.016 0.022 0.028
## mode kld
## (Intercept) -0.291 0
## vandix_z 0.185 0
## ln_roads_km_c 0.844 0
## roads_prop_highway_arterial_z 0.301 0
## ale_index_z 0.041 0
## canbics_index_z 0.092 0
## population_100_c 0.022 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.822 0.069 0.692 0.819 0.964 0.817
## Phi for ui 0.871 0.036 0.790 0.874 0.931 0.880
##
## Deviance Information Criterion (DIC) ...............: 8930.87
## Deviance Information Criterion (DIC, saturated) ....: 4377.08
## Effective number of parameters .....................: 1046.94
##
## Watanabe-Akaike information criterion (WAIC) ...: 8897.21
## Effective number of parameters .................: 790.56
##
## Marginal log-Likelihood: -3706.06
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 3: Total Casualty Cyclist Claim Crashes
van_models_3 <- cma_models(vancouver_sp@data,"Vancouver-Squamish","n_pedestrian_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"van_pedestrian <- clean_fixed_effects(van_models_3$fixed_effects)
map(van_models_3$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 18, Running = 4.19, Post = 0.429, Total = 22.7
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -0.035 0.027 -0.089 -0.034 0.018 -0.034 0
## vandix_z 0.262 0.036 0.191 0.262 0.333 0.262 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.496 0.035 0.430 0.494 0.567 0.492
## Phi for ui 0.663 0.044 0.573 0.664 0.746 0.666
##
## Deviance Information Criterion (DIC) ...............: 11940.16
## Deviance Information Criterion (DIC, saturated) ....: 5994.58
## Effective number of parameters .....................: 2208.49
##
## Watanabe-Akaike information criterion (WAIC) ...: 12274.22
## Effective number of parameters .................: 1755.74
##
## Marginal log-Likelihood: -5474.94
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 18.5, Running = 4.56, Post = 0.393, Total = 23.5
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 0.319 0.022 0.275 0.319 0.363 0.319 0
## vandix_z 0.304 0.032 0.242 0.304 0.367 0.304 0
## ln_roads_km_c 0.997 0.034 0.930 0.997 1.065 0.997 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.45 0.028 0.396 0.450 0.506 0.450
## Phi for ui 0.95 0.019 0.907 0.953 0.979 0.957
##
## Deviance Information Criterion (DIC) ...............: 11175.53
## Deviance Information Criterion (DIC, saturated) ....: 5229.93
## Effective number of parameters .....................: 1630.40
##
## Watanabe-Akaike information criterion (WAIC) ...: 11151.04
## Effective number of parameters .................: 1185.45
##
## Marginal log-Likelihood: -5141.60
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 18.5, Running = 5.45, Post = 0.451, Total = 24.4
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode
## (Intercept) 0.115 0.028 0.060 0.115 0.168 0.115
## vandix_z 0.268 0.031 0.208 0.268 0.328 0.268
## ln_roads_km_c 0.717 0.039 0.641 0.717 0.793 0.717
## roads_prop_highway_arterial_z 0.394 0.022 0.351 0.394 0.438 0.394
## ale_index_z 0.232 0.046 0.141 0.232 0.322 0.232
## canbics_index_z 0.042 0.047 -0.051 0.042 0.135 0.042
## population_100_c 0.027 0.003 0.021 0.027 0.034 0.027
## kld
## (Intercept) 0
## vandix_z 0
## ln_roads_km_c 0
## roads_prop_highway_arterial_z 0
## ale_index_z 0
## canbics_index_z 0
## population_100_c 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.713 0.059 0.602 0.712 0.835 0.709
## Phi for ui 0.850 0.039 0.763 0.853 0.917 0.860
##
## Deviance Information Criterion (DIC) ...............: 11047.95
## Deviance Information Criterion (DIC, saturated) ....: 5102.36
## Effective number of parameters .....................: 1480.57
##
## Watanabe-Akaike information criterion (WAIC) ...: 10986.17
## Effective number of parameters .................: 1069.75
##
## Marginal log-Likelihood: -4953.61
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')van_results <- bind_rows(van_all %>% mutate(Region = "Vancouver-Squamish",Outcome = "All Injury Claims"),
van_cyclist %>% mutate(Region = "Vancouver-Squamish",Outcome = "Cyclist Injury Claims"),
van_pedestrian %>% mutate(Region = "Vancouver-Squamish",Outcome = "Pedestrian Injury Claims")
) %>%
filter(variable == "vandix_z") %>%
select(Region,Outcome,everything())Victoria
Descriptive Statistics
region | variable | n | sum | mean | sd | min | max | missing | n_zero | p_zero |
|---|---|---|---|---|---|---|---|---|---|---|
Victoria | n_claims | 581 | 67,182.0 | 115.6 | 176.1 | 0.0 | 1,665.0 | 0 | 3 | 0.5 |
n_casualty_claims | 12,338.0 | 21.2 | 34.2 | 0.0 | 305.0 | 0 | 26 | 4.5 | ||
n_cyclist_claims | 1,471.0 | 2.5 | 4.4 | 0.0 | 53.0 | 0 | 196 | 33.7 | ||
n_cyclist_casualty_claims | 1,027.0 | 1.8 | 3.2 | 0.0 | 33.0 | 0 | 253 | 43.5 | ||
n_pedestrian_claims | 896.0 | 1.5 | 3.4 | 0.0 | 49.0 | 0 | 282 | 48.5 | ||
n_pedestrian_casualty_claims | 715.0 | 1.2 | 2.5 | 0.0 | 37.0 | 0 | 304 | 52.3 | ||
population | 397,237.0 | 683.7 | 432.5 | 0.0 | 4,739.0 | 0 | 3 | 0.5 | ||
total_roads_km | 3,261.2 | 5.6 | 5.6 | 0.0 | 48.0 | 0 | 2 | 0.3 | ||
roads_prop_highway_arterial | 67.1 | 0.1 | 0.1 | 0.0 | 0.7 | 2 | 238 | 41.0 | ||
no_highschool_prevalance | 66.7 | 0.1 | 0.1 | 0.0 | 0.4 | 30 | 0 | 0.0 | ||
unemployment_rate | 3,099.4 | 5.6 | 3.9 | 0.0 | 30.0 | 30 | 78 | 13.4 | ||
hh_avg_income | 41,608,895.0 | 77,340.0 | 26,930.8 | 26,088.0 | 233,035.0 | 43 | 0 | 0.0 | ||
participation_rate | 35,068.6 | 63.6 | 11.1 | 18.6 | 92.7 | 30 | 0 | 0.0 | ||
university_degree_prevalance | 329.1 | 0.6 | 0.1 | 0.3 | 0.9 | 30 | 0 | 0.0 | ||
lone_parent_fam_prevalence | 84.0 | 0.2 | 0.1 | 0.0 | 0.6 | 32 | 5 | 0.9 | ||
home_owner_prevalence | 360.4 | 0.7 | 0.2 | 0.0 | 1.0 | 30 | 3 | 0.5 | ||
vandix | -130.2 | -0.2 | 0.5 | -1.5 | 3.0 | 45 | 0 | 0.0 | ||
ale_index | 335.8 | 0.6 | 1.8 | -2.1 | 7.1 | 28 | 2 | 0.3 | ||
canbics_index | 1,170.3 | 2.1 | 1.7 | 0.0 | 8.2 | 28 | 69 | 11.9 |
Maps
claims <- ggplot() +
geom_sf(data = victoria_sf, aes(fill = n_casualty_claims,colour=n_casualty_claims)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
theme_void() +
ggtitle(
"Vancouver"
)
vandix <- ggplot() +
geom_sf(data = victoria_sf, aes(fill = vandix_z_c,colour=vandix_z_c)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
scale_colour_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
theme_void()
total_roads <- ggplot() +
geom_sf(data = victoria_sf, aes(fill = total_roads_km,colour=total_roads_km)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
theme_void()
cowplot::plot_grid(claims,vandix,total_roads,ncol=1)
Setting up Spatial Weights Matrix
#### Define Spatial Neighrbourhoods
victoria_sp <- as(victoria_sf, "Spatial")
victoria_sp$ui <- 1:nrow(victoria_sp@data)
coords <- coordinates(victoria_sp)
victoria_nb <- poly2nb(victoria_sp, queen = TRUE)
victoria_nb## Neighbour list object:
## Number of regions: 581
## Number of nonzero links: 3468
## Percentage nonzero weights: 1.02737
## Average number of links: 5.969019#assign nearest neighbour for no links
# victoria_nb <- assign_nearest_neighbors(victoria_nb,victoria_sp)
plot(victoria_sp, border = grey(0.5))
plot(victoria_nb,
coords = coords,
add = TRUE, pch = 16, lwd = 2)
Spatial Dependency in each crash type
listw <- nb2listw(victoria_nb,zero.policy = TRUE)
all_cc_mi <- moran.test(victoria_sf$n_casualty_claims, listw,zero.policy = TRUE)
all_cc_mi##
## Moran I test under randomisation
##
## data: victoria_sf$n_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 16.303, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.3913583171 -0.0017241379 0.0005813523cyc_cc_mi <- moran.test(victoria_sf$n_cyclist_casualty_claims, listw,zero.policy = TRUE)
cyc_cc_mi##
## Moran I test under randomisation
##
## data: victoria_sf$n_cyclist_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 13.667, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.3259641289 -0.0017241379 0.0005749071pd_cc_mi <- moran.test(victoria_sf$n_pedestrian_casualty_claims, listw,zero.policy = TRUE)
pd_cc_mi##
## Moran I test under randomisation
##
## data: victoria_sf$n_pedestrian_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 15.441, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.3531561826 -0.0017241379 0.0005282317BYM2 Models
nb2INLA("victoria.adj", victoria_nb)
g <- inla.read.graph(filename = "victoria.adj")
# Model Set 1: Total Casualty Claim Crashes
victoria_models_1 <- cma_models(victoria_sp@data, "Victoria", "n_casualty_claims", "vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"victoria_all <- clean_fixed_effects(victoria_models_1$fixed_effects)
map(victoria_models_1$models, ~ summary(.x))## $Unadjusted
## Time used:
## Pre = 14.2, Running = 0.866, Post = 0.148, Total = 15.3
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 2.339 0.035 2.270 2.339 2.409 2.339 0
## vandix_z 0.196 0.065 0.068 0.196 0.325 0.196 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.490 0.053 0.390 0.488 0.599 0.486
## Phi for ui 0.857 0.059 0.721 0.865 0.947 0.882
##
## Deviance Information Criterion (DIC) ...............: 3477.53
## Deviance Information Criterion (DIC, saturated) ....: 1131.62
## Effective number of parameters .....................: 510.12
##
## Watanabe-Akaike information criterion (WAIC) ...: 3449.37
## Effective number of parameters .................: 336.48
##
## Marginal log-Likelihood: -1963.43
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.2, Running = 0.894, Post = 0.162, Total = 15.2
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 2.361 0.025 2.311 2.361 2.41 2.361 0
## vandix_z 0.326 0.053 0.222 0.326 0.43 0.326 0
## ln_roads_km_c 1.237 0.073 1.095 1.237 1.38 1.237 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.642 0.050 0.547 0.640 0.744 0.639
## Phi for ui 0.981 0.019 0.930 0.987 0.999 0.996
##
## Deviance Information Criterion (DIC) ...............: 3406.34
## Deviance Information Criterion (DIC, saturated) ....: 1060.42
## Effective number of parameters .....................: 465.82
##
## Watanabe-Akaike information criterion (WAIC) ...: 3339.80
## Effective number of parameters .................: 285.81
##
## Marginal log-Likelihood: -1849.39
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.1, Running = 0.947, Post = 0.163, Total = 15.2
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode
## (Intercept) 2.514 0.046 2.424 2.515 2.604 2.515
## vandix_z 0.202 0.050 0.104 0.202 0.301 0.202
## ln_roads_km_c 1.007 0.080 0.850 1.007 1.165 1.007
## roads_prop_highway_arterial_z 0.462 0.049 0.367 0.462 0.558 0.462
## ale_index_z 0.409 0.168 0.081 0.409 0.739 0.409
## canbics_index_z 0.340 0.156 0.032 0.340 0.646 0.340
## population_100_c 0.015 0.009 -0.002 0.015 0.031 0.015
## kld
## (Intercept) 0
## vandix_z 0
## ln_roads_km_c 0
## roads_prop_highway_arterial_z 0
## ale_index_z 0
## canbics_index_z 0
## population_100_c 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.845 0.087 0.679 0.844 1.021 0.847
## Phi for ui 0.948 0.042 0.838 0.959 0.994 0.982
##
## Deviance Information Criterion (DIC) ...............: 3403.57
## Deviance Information Criterion (DIC, saturated) ....: 1057.64
## Effective number of parameters .....................: 453.65
##
## Watanabe-Akaike information criterion (WAIC) ...: 3344.18
## Effective number of parameters .................: 283.52
##
## Marginal log-Likelihood: -1822.76
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 2: Total Casualty Cyclist Claim Crashes
victoria_models_2 <- cma_models(victoria_sp@data,"Victoria","n_cyclist_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"victoria_cyclist <- clean_fixed_effects(victoria_models_2$fixed_effects)
map(victoria_models_2$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 13.9, Running = 0.818, Post = 0.145, Total = 14.8
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -0.141 0.070 -0.282 -0.14 -0.006 -0.14 0
## vandix_z 0.150 0.085 -0.017 0.15 0.318 0.15 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.679 0.094 0.511 0.673 0.882 0.661
## Phi for ui 0.623 0.101 0.414 0.628 0.806 0.637
##
## Deviance Information Criterion (DIC) ...............: 1756.63
## Deviance Information Criterion (DIC, saturated) ....: 872.31
## Effective number of parameters .....................: 289.24
##
## Watanabe-Akaike information criterion (WAIC) ...: 1758.23
## Effective number of parameters .................: 211.24
##
## Marginal log-Likelihood: -721.47
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14, Running = 0.803, Post = 0.157, Total = 15
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -0.144 0.065 -0.273 -0.143 -0.018 -0.143 0
## vandix_z 0.250 0.080 0.094 0.250 0.406 0.250 0
## ln_roads_km_c 1.096 0.115 0.871 1.095 1.322 1.095 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.663 0.096 0.489 0.658 0.868 0.650
## Phi for ui 0.902 0.064 0.739 0.917 0.984 0.951
##
## Deviance Information Criterion (DIC) ...............: 1693.39
## Deviance Information Criterion (DIC, saturated) ....: 809.07
## Effective number of parameters .....................: 240.45
##
## Watanabe-Akaike information criterion (WAIC) ...: 1684.94
## Effective number of parameters .................: 174.53
##
## Marginal log-Likelihood: -684.76
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.2, Running = 0.945, Post = 0.152, Total = 15.3
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode
## (Intercept) 0.039 0.080 -0.120 0.040 0.193 0.040
## vandix_z 0.118 0.080 -0.039 0.118 0.275 0.118
## ln_roads_km_c 0.887 0.129 0.635 0.887 1.140 0.887
## roads_prop_highway_arterial_z 0.308 0.073 0.165 0.308 0.451 0.308
## ale_index_z 0.401 0.222 -0.033 0.401 0.838 0.401
## canbics_index_z 0.587 0.213 0.169 0.588 1.006 0.588
## population_100_c 0.031 0.014 0.004 0.031 0.058 0.031
## kld
## (Intercept) 0
## vandix_z 0
## ln_roads_km_c 0
## roads_prop_highway_arterial_z 0
## ale_index_z 0
## canbics_index_z 0
## population_100_c 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.867 0.138 0.622 0.858 1.164 0.842
## Phi for ui 0.828 0.089 0.616 0.843 0.957 0.880
##
## Deviance Information Criterion (DIC) ...............: 1686.04
## Deviance Information Criterion (DIC, saturated) ....: 801.72
## Effective number of parameters .....................: 225.96
##
## Watanabe-Akaike information criterion (WAIC) ...: 1679.77
## Effective number of parameters .................: 167.01
##
## Marginal log-Likelihood: -686.41
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 3: Total Casualty Cyclist Claim Crashes
victoria_models_3 <- cma_models(victoria_sp@data,"Victoria","n_pedestrian_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"victoria_pedestrian <- clean_fixed_effects(victoria_models_3$fixed_effects)
map(victoria_models_3$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 14.1, Running = 0.812, Post = 0.145, Total = 15.1
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -0.439 0.074 -0.589 -0.438 -0.296 -0.439 0
## vandix_z 0.350 0.084 0.184 0.349 0.516 0.349 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.633 0.099 0.458 0.626 0.847 0.615
## Phi for ui 0.791 0.096 0.569 0.805 0.937 0.837
##
## Deviance Information Criterion (DIC) ...............: 1480.09
## Deviance Information Criterion (DIC, saturated) ....: 770.60
## Effective number of parameters .....................: 228.04
##
## Watanabe-Akaike information criterion (WAIC) ...: 1476.73
## Effective number of parameters .................: 167.60
##
## Marginal log-Likelihood: -560.34
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.2, Running = 0.812, Post = 0.15, Total = 15.2
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -0.430 0.071 -0.573 -0.429 -0.292 -0.429 0
## vandix_z 0.452 0.081 0.292 0.452 0.612 0.452 0
## ln_roads_km_c 0.854 0.118 0.625 0.854 1.086 0.854 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.619 0.091 0.458 0.613 0.816 0.602
## Phi for ui 0.950 0.047 0.822 0.964 0.996 0.988
##
## Deviance Information Criterion (DIC) ...............: 1438.88
## Deviance Information Criterion (DIC, saturated) ....: 729.38
## Effective number of parameters .....................: 198.78
##
## Watanabe-Akaike information criterion (WAIC) ...: 1428.78
## Effective number of parameters .................: 144.38
##
## Marginal log-Likelihood: -540.86
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.2, Running = 1.56, Post = 1.14, Total = 16.9
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -0.313 0.090 -0.493 -0.312 -0.141
## vandix_z 0.317 0.082 0.156 0.317 0.478
## ln_roads_km_c 0.578 0.135 0.313 0.578 0.845
## roads_prop_highway_arterial_z 0.384 0.079 0.228 0.383 0.539
## ale_index_z 0.540 0.237 0.080 0.539 1.008
## canbics_index_z 0.326 0.235 -0.136 0.327 0.786
## population_100_c 0.037 0.014 0.010 0.037 0.064
## mode kld
## (Intercept) -0.312 0
## vandix_z 0.317 0
## ln_roads_km_c 0.578 0
## roads_prop_highway_arterial_z 0.383 0
## ale_index_z 0.539 0
## canbics_index_z 0.327 0
## population_100_c 0.037 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.758 0.122 0.543 0.750 1.023 0.734
## Phi for ui 0.924 0.066 0.745 0.944 0.993 0.980
##
## Deviance Information Criterion (DIC) ...............: 1426.26
## Deviance Information Criterion (DIC, saturated) ....: 716.76
## Effective number of parameters .....................: 185.00
##
## Watanabe-Akaike information criterion (WAIC) ...: 1417.96
## Effective number of parameters .................: 136.61
##
## Marginal log-Likelihood: -542.77
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')victoria_results <- bind_rows(victoria_all %>% mutate(Region = "Victoria",Outcome = "All Injury Claims"),
victoria_cyclist %>% mutate(Region = "Victoria",Outcome = "Cyclist Injury Claims"),
victoria_pedestrian %>% mutate(Region = "Victoria",Outcome = "Pedestrian Injury Claims")
) %>%
filter(variable == "vandix_z") %>%
select(Region,Outcome,everything())Central Island - Powell River
Descriptive Statistics
region | variable | n | sum | mean | sd | min | max | missing | n_zero | p_zero |
|---|---|---|---|---|---|---|---|---|---|---|
Central Island-Powell River | n_claims | 611 | 62,992.0 | 103.1 | 164.7 | 0.0 | 1,934.0 | 0 | 17 | 2.8 |
n_casualty_claims | 10,762.0 | 17.6 | 28.3 | 0.0 | 203.0 | 0 | 47 | 7.7 | ||
n_cyclist_claims | 475.0 | 0.8 | 1.5 | 0.0 | 14.0 | 0 | 365 | 59.7 | ||
n_cyclist_casualty_claims | 271.0 | 0.4 | 0.9 | 0.0 | 10.0 | 0 | 435 | 71.2 | ||
n_pedestrian_claims | 673.0 | 1.1 | 2.4 | 0.0 | 23.0 | 0 | 361 | 59.1 | ||
n_pedestrian_casualty_claims | 547.0 | 0.9 | 2.0 | 0.0 | 18.0 | 0 | 387 | 63.3 | ||
population | 357,193.0 | 584.6 | 307.8 | 0.0 | 2,466.0 | 0 | 14 | 2.3 | ||
total_roads_km | 5,652.1 | 9.3 | 9.9 | 0.0 | 86.5 | 0 | 10 | 1.6 | ||
roads_prop_highway_arterial | 78.2 | 0.1 | 0.2 | 0.0 | 0.8 | 10 | 262 | 42.9 | ||
no_highschool_prevalance | 101.2 | 0.2 | 0.1 | 0.0 | 0.6 | 38 | 1 | 0.2 | ||
unemployment_rate | 4,784.2 | 8.3 | 5.2 | 0.0 | 40.0 | 38 | 43 | 7.0 | ||
hh_avg_income | 36,682,051.0 | 65,738.4 | 17,197.1 | 29,925.0 | 211,531.0 | 53 | 0 | 0.0 | ||
participation_rate | 32,393.0 | 56.5 | 11.0 | 12.0 | 86.4 | 38 | 0 | 0.0 | ||
university_degree_prevalance | 300.2 | 0.5 | 0.1 | 0.2 | 0.8 | 38 | 0 | 0.0 | ||
lone_parent_fam_prevalence | 95.7 | 0.2 | 0.1 | 0.0 | 0.6 | 38 | 3 | 0.5 | ||
home_owner_prevalence | 426.2 | 0.7 | 0.2 | 0.0 | 1.0 | 38 | 1 | 0.2 | ||
vandix | 111.8 | 0.2 | 0.7 | -1.1 | 4.0 | 53 | 0 | 0.0 | ||
ale_index | -590.4 | -1.1 | 0.7 | -2.1 | 1.5 | 57 | 0 | 0.0 | ||
canbics_index | 326.8 | 0.6 | 0.9 | 0.0 | 4.2 | 55 | 241 | 39.4 |
Maps
claims <- ggplot() +
geom_sf(data = nanaimo_sf, aes(fill = n_casualty_claims,colour=n_casualty_claims)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
theme_void() +
ggtitle(
"Vancouver"
)
vandix <- ggplot() +
geom_sf(data = nanaimo_sf, aes(fill = vandix_z_c,colour=vandix_z_c)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
scale_colour_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
theme_void()
total_roads <- ggplot() +
geom_sf(data = nanaimo_sf, aes(fill = total_roads_km,colour=total_roads_km)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
theme_void()
cowplot::plot_grid(claims,vandix,total_roads,ncol=1)
Setting up Spatial Weights Matrix
#### Define Spatial Neighrbourhoods
nanaimo_sp <- as(nanaimo_sf, "Spatial")
nanaimo_sp$ui <- 1:nrow(nanaimo_sp@data)
coords <- coordinates(nanaimo_sp)
nanaimo_nb <- poly2nb(nanaimo_sp, queen = TRUE)## Warning in poly2nb(nanaimo_sp, queen = TRUE): some observations have no neighbours;
## if this seems unexpected, try increasing the snap argument.## Warning in poly2nb(nanaimo_sp, queen = TRUE): neighbour object has 6 sub-graphs;
## if this sub-graph count seems unexpected, try increasing the snap argument.nanaimo_nb## Neighbour list object:
## Number of regions: 611
## Number of nonzero links: 3330
## Percentage nonzero weights: 0.8919938
## Average number of links: 5.450082
## 1 region with no links:
## 169
## 6 disjoint connected subgraphs#assign nearest neighbour for no links
nanaimo_nb <- assign_nearest_neighbors(nanaimo_nb,nanaimo_sp)
plot(nanaimo_sp, border = grey(0.5))
plot(nanaimo_nb,
coords = coords,
add = TRUE, pch = 16, lwd = 2)
Spatial Dependency in each crash type
listw <- nb2listw(nanaimo_nb,zero.policy = TRUE)
all_cc_mi <- moran.test(nanaimo_sf$n_casualty_claims, listw,zero.policy = TRUE)
all_cc_mi##
## Moran I test under randomisation
##
## data: nanaimo_sf$n_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 12.06, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.3077074730 -0.0016393443 0.0006579567cyc_cc_mi <- moran.test(nanaimo_sf$n_cyclist_casualty_claims, listw,zero.policy = TRUE)
cyc_cc_mi##
## Moran I test under randomisation
##
## data: nanaimo_sf$n_cyclist_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 6.5732, p-value = 2.463e-11
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.1652016950 -0.0016393443 0.0006442558pd_cc_mi <- moran.test(nanaimo_sf$n_pedestrian_casualty_claims, listw,zero.policy = TRUE)
pd_cc_mi##
## Moran I test under randomisation
##
## data: nanaimo_sf$n_pedestrian_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 8.4624, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.2136910003 -0.0016393443 0.0006474725BYM2 Models
nb2INLA("nanaimo.adj", nanaimo_nb)
g <- inla.read.graph(filename = "nanaimo.adj")
# Model Set 1: Total Casualty Claim Crashes
nanaimo_models_1 <- cma_models(nanaimo_sp@data, "Central Island-Powell River", "n_casualty_claims", "vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"nanaimo_all <- clean_fixed_effects(nanaimo_models_1$fixed_effects)
map(nanaimo_models_1$models, ~ summary(.x))## $Unadjusted
## Time used:
## Pre = 15.4, Running = 0.911, Post = 0.19, Total = 16.5
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 1.926 0.050 1.826 1.926 2.024 1.926 0
## vandix_z 0.295 0.061 0.175 0.295 0.416 0.295 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.410 0.058 0.309 0.406 0.537 0.397
## Phi for ui 0.571 0.092 0.383 0.574 0.741 0.583
##
## Deviance Information Criterion (DIC) ...............: 3536.92
## Deviance Information Criterion (DIC, saturated) ....: 1244.46
## Effective number of parameters .....................: 555.03
##
## Watanabe-Akaike information criterion (WAIC) ...: 3574.86
## Effective number of parameters .................: 405.59
##
## Marginal log-Likelihood: -2193.01
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 15.9, Running = 1.18, Post = 0.168, Total = 17.2
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 1.357 0.060 1.238 1.357 1.475 1.357 0
## vandix_z 0.435 0.053 0.330 0.435 0.539 0.435 0
## ln_roads_km_c 1.024 0.069 0.888 1.023 1.160 1.023 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.442 0.054 0.345 0.438 0.558 0.432
## Phi for ui 0.736 0.057 0.613 0.740 0.836 0.748
##
## Deviance Information Criterion (DIC) ...............: 3467.76
## Deviance Information Criterion (DIC, saturated) ....: 1175.29
## Effective number of parameters .....................: 523.49
##
## Watanabe-Akaike information criterion (WAIC) ...: 3456.68
## Effective number of parameters .................: 356.15
##
## Marginal log-Likelihood: -2098.24
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 15.6, Running = 1.08, Post = 0.193, Total = 16.9
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) 1.961 0.117 1.730 1.962 2.190
## vandix_z 0.317 0.046 0.227 0.317 0.408
## ln_roads_km_c 0.650 0.073 0.506 0.650 0.794
## roads_prop_highway_arterial_z 0.668 0.047 0.575 0.668 0.761
## ale_index_z 1.419 0.246 0.937 1.419 1.903
## canbics_index_z -0.780 0.231 -1.233 -0.779 -0.328
## population_100_c 0.063 0.013 0.037 0.063 0.090
## mode kld
## (Intercept) 1.962 0
## vandix_z 0.317 0
## ln_roads_km_c 0.650 0
## roads_prop_highway_arterial_z 0.668 0
## ale_index_z 1.419 0
## canbics_index_z -0.779 0
## population_100_c 0.063 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.670 0.085 0.517 0.665 0.852 0.656
## Phi for ui 0.726 0.060 0.598 0.730 0.832 0.738
##
## Deviance Information Criterion (DIC) ...............: 3436.59
## Deviance Information Criterion (DIC, saturated) ....: 1144.12
## Effective number of parameters .....................: 492.67
##
## Watanabe-Akaike information criterion (WAIC) ...: 3402.21
## Effective number of parameters .................: 324.32
##
## Marginal log-Likelihood: -2024.42
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 2: Total Casualty Cyclist Claim Crashes
nanaimo_models_2 <- cma_models(nanaimo_sp@data,"Central Island-Powell River","n_cyclist_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"nanaimo_cyclist <- clean_fixed_effects(nanaimo_models_2$fixed_effects)
map(nanaimo_models_2$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 15.6, Running = 0.969, Post = 0.154, Total = 16.8
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.429 0.110 -1.655 -1.427 -1.220 -1.427 0
## vandix_z 0.335 0.075 0.188 0.335 0.483 0.335 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1.06 0.236 0.680 1.031 1.605 0.972
## Phi for ui 0.32 0.142 0.089 0.305 0.626 0.266
##
## Deviance Information Criterion (DIC) ...............: 1012.10
## Deviance Information Criterion (DIC, saturated) ....: 614.45
## Effective number of parameters .....................: 137.56
##
## Watanabe-Akaike information criterion (WAIC) ...: 1008.51
## Effective number of parameters .................: 106.42
##
## Marginal log-Likelihood: -435.69
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 15.3, Running = 1.01, Post = 0.158, Total = 16.5
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.731 0.139 -2.012 -1.729 -1.465 -1.729 0
## vandix_z 0.399 0.076 0.250 0.399 0.548 0.399 0
## ln_roads_km_c 0.506 0.117 0.280 0.505 0.737 0.505 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.927 0.231 0.561 0.898 1.463 0.841
## Phi for ui 0.551 0.151 0.249 0.557 0.823 0.577
##
## Deviance Information Criterion (DIC) ...............: 993.69
## Deviance Information Criterion (DIC, saturated) ....: 596.04
## Effective number of parameters .....................: 126.27
##
## Watanabe-Akaike information criterion (WAIC) ...: 993.42
## Effective number of parameters .................: 100.60
##
## Marginal log-Likelihood: -431.79
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 15.2, Running = 1.09, Post = 0.173, Total = 16.5
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -1.152 0.208 -1.568 -1.149 -0.753
## vandix_z 0.354 0.077 0.203 0.354 0.504
## ln_roads_km_c 0.292 0.140 0.018 0.291 0.568
## roads_prop_highway_arterial_z 0.431 0.082 0.272 0.431 0.592
## ale_index_z 1.156 0.413 0.350 1.154 1.969
## canbics_index_z -0.321 0.366 -1.042 -0.319 0.394
## population_100_c 0.060 0.024 0.014 0.060 0.107
## mode kld
## (Intercept) -1.149 0
## vandix_z 0.354 0
## ln_roads_km_c 0.291 0
## roads_prop_highway_arterial_z 0.431 0
## ale_index_z 1.154 0
## canbics_index_z -0.319 0
## population_100_c 0.060 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1.177 0.313 0.691 1.133 1.915 1.05
## Phi for ui 0.496 0.156 0.197 0.497 0.788 0.51
##
## Deviance Information Criterion (DIC) ...............: 971.41
## Deviance Information Criterion (DIC, saturated) ....: 573.76
## Effective number of parameters .....................: 113.46
##
## Watanabe-Akaike information criterion (WAIC) ...: 972.04
## Effective number of parameters .................: 92.21
##
## Marginal log-Likelihood: -433.56
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 3: Total Casualty Cyclist Claim Crashes
nanaimo_models_3 <- cma_models(nanaimo_sp@data,"Central Island-Powell River","n_pedestrian_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"nanaimo_pedestrian <- clean_fixed_effects(nanaimo_models_3$fixed_effects)
map(nanaimo_models_3$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 15.2, Running = 1.02, Post = 0.159, Total = 16.3
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.252 0.116 -1.494 -1.248 -1.036 -1.248 0
## vandix_z 0.556 0.078 0.405 0.556 0.710 0.556 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.605 0.100 0.439 0.595 0.830 0.572
## Phi for ui 0.206 0.128 0.030 0.180 0.511 0.102
##
## Deviance Information Criterion (DIC) ...............: 1356.86
## Deviance Information Criterion (DIC, saturated) ....: 791.04
## Effective number of parameters .....................: 274.96
##
## Watanabe-Akaike information criterion (WAIC) ...: 1378.25
## Effective number of parameters .................: 208.18
##
## Marginal log-Likelihood: -648.63
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 15.5, Running = 1.05, Post = 0.162, Total = 16.8
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.598 0.140 -1.882 -1.595 -1.330 -1.595 0
## vandix_z 0.611 0.080 0.456 0.611 0.769 0.611 0
## ln_roads_km_c 0.653 0.129 0.401 0.652 0.908 0.652 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.468 0.102 0.302 0.456 0.701 0.433
## Phi for ui 0.587 0.140 0.295 0.596 0.830 0.627
##
## Deviance Information Criterion (DIC) ...............: 1314.98
## Deviance Information Criterion (DIC, saturated) ....: 749.16
## Effective number of parameters .....................: 241.00
##
## Watanabe-Akaike information criterion (WAIC) ...: 1327.06
## Effective number of parameters .................: 182.32
##
## Marginal log-Likelihood: -641.60
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 15.4, Running = 1.04, Post = 0.155, Total = 16.5
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -0.755 0.205 -1.167 -0.752 -0.361
## vandix_z 0.518 0.075 0.370 0.518 0.666
## ln_roads_km_c 0.328 0.141 0.052 0.327 0.605
## roads_prop_highway_arterial_z 0.552 0.081 0.394 0.551 0.713
## ale_index_z 2.371 0.393 1.605 2.369 3.147
## canbics_index_z -1.147 0.349 -1.838 -1.146 -0.468
## population_100_c 0.115 0.023 0.071 0.115 0.160
## mode kld
## (Intercept) -0.752 0
## vandix_z 0.518 0
## ln_roads_km_c 0.327 0
## roads_prop_highway_arterial_z 0.551 0
## ale_index_z 2.369 0
## canbics_index_z -1.145 0
## population_100_c 0.115 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.749 0.161 0.489 0.731 1.119 0.692
## Phi for ui 0.446 0.150 0.166 0.444 0.736 0.451
##
## Deviance Information Criterion (DIC) ...............: 1265.20
## Deviance Information Criterion (DIC, saturated) ....: 699.38
## Effective number of parameters .....................: 201.41
##
## Watanabe-Akaike information criterion (WAIC) ...: 1265.47
## Effective number of parameters .................: 150.79
##
## Marginal log-Likelihood: -617.89
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')nanaimo_results <- bind_rows(nanaimo_all %>% mutate(Region = "Central Island-Powell River",Outcome = "All Injury Claims"),
nanaimo_cyclist %>% mutate(Region = "Central Island-Powell River",Outcome = "Cyclist Injury Claims"),
nanaimo_pedestrian %>% mutate(Region = "Central Island-Powell River",Outcome = "Pedestrian Injury Claims")
) %>%
filter(variable == "vandix_z") %>%
select(Region,Outcome,everything())Okanagan
Descriptive Statistics
region | variable | n | sum | mean | sd | min | max | missing | n_zero | p_zero |
|---|---|---|---|---|---|---|---|---|---|---|
Okanagan | n_claims | 466 | 70,024.0 | 150.3 | 325.5 | 0.0 | 4,324.0 | 0 | 10 | 2.1 |
n_casualty_claims | 13,236.0 | 28.4 | 66.6 | 0.0 | 938.0 | 0 | 38 | 8.2 | ||
n_cyclist_claims | 761.0 | 1.6 | 3.8 | 0.0 | 38.0 | 0 | 252 | 54.1 | ||
n_cyclist_casualty_claims | 485.0 | 1.0 | 2.3 | 0.0 | 23.0 | 0 | 286 | 61.4 | ||
n_pedestrian_claims | 696.0 | 1.5 | 3.9 | 0.0 | 41.0 | 0 | 259 | 55.6 | ||
n_pedestrian_casualty_claims | 564.0 | 1.2 | 3.2 | 0.0 | 35.0 | 0 | 274 | 58.8 | ||
population | 336,628.0 | 722.4 | 489.5 | 0.0 | 4,622.0 | 0 | 10 | 2.1 | ||
total_roads_km | 4,759.1 | 10.2 | 13.6 | 0.0 | 139.3 | 0 | 2 | 0.4 | ||
roads_prop_highway_arterial | 60.1 | 0.1 | 0.2 | 0.0 | 1.0 | 2 | 199 | 42.7 | ||
no_highschool_prevalance | 62.0 | 0.2 | 0.1 | 0.0 | 0.4 | 88 | 4 | 0.9 | ||
unemployment_rate | 2,933.5 | 7.8 | 4.7 | 0.0 | 33.3 | 88 | 28 | 6.0 | ||
hh_avg_income | 26,784,003.0 | 72,980.9 | 26,808.0 | 26,241.0 | 247,493.0 | 99 | 0 | 0.0 | ||
participation_rate | 23,166.9 | 61.3 | 12.9 | 9.2 | 91.5 | 88 | 0 | 0.0 | ||
university_degree_prevalance | 199.0 | 0.5 | 0.1 | 0.3 | 0.8 | 88 | 0 | 0.0 | ||
lone_parent_fam_prevalence | 59.9 | 0.2 | 0.1 | 0.0 | 0.5 | 88 | 4 | 0.9 | ||
home_owner_prevalence | 275.1 | 0.7 | 0.2 | 0.0 | 1.3 | 88 | 0 | 0.0 | ||
vandix | 42.6 | 0.1 | 0.6 | -1.2 | 2.2 | 99 | 0 | 0.0 | ||
ale_index | -314.5 | -0.8 | 1.0 | -2.1 | 2.3 | 93 | 1 | 0.2 | ||
canbics_index | 709.8 | 1.9 | 1.8 | 0.0 | 7.2 | 93 | 82 | 17.6 |
Maps
claims <- ggplot() +
geom_sf(data = okanagan_sf, aes(fill = n_casualty_claims,colour=n_casualty_claims)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
theme_void() +
ggtitle(
"Vancouver"
)
vandix <- ggplot() +
geom_sf(data = okanagan_sf, aes(fill = vandix_z_c,colour=vandix_z_c)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
scale_colour_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
theme_void()
total_roads <- ggplot() +
geom_sf(data = okanagan_sf, aes(fill = total_roads_km,colour=total_roads_km)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
theme_void()
cowplot::plot_grid(claims,vandix,total_roads,ncol=1)
Setting up Spatial Weights Matrix
#### Define Spatial Neighrbourhoods
okanagan_sp <- as(okanagan_sf, "Spatial")
okanagan_sp$ui <- 1:nrow(okanagan_sp@data)
coords <- coordinates(okanagan_sp)
okanagan_nb <- poly2nb(okanagan_sp, queen = TRUE,snap = 1)
okanagan_nb## Neighbour list object:
## Number of regions: 466
## Number of nonzero links: 2728
## Percentage nonzero weights: 1.25624
## Average number of links: 5.854077#assign nearest neighbour for no links
# okanagan_nb <- assign_nearest_neighbors(okanagan_nb,okanagan_sp)
plot(okanagan_sp, border = grey(0.5))
plot(okanagan_nb,
coords = coords,
add = TRUE, pch = 16, lwd = 2)
Spatial Dependency in each crash type
listw <- nb2listw(okanagan_nb,zero.policy = TRUE)
all_cc_mi <- moran.test(okanagan_sf$n_casualty_claims, listw,zero.policy = TRUE)
all_cc_mi##
## Moran I test under randomisation
##
## data: okanagan_sf$n_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 14.429, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.3626842832 -0.0021505376 0.0006393647cyc_cc_mi <- moran.test(okanagan_sf$n_cyclist_casualty_claims, listw,zero.policy = TRUE)
cyc_cc_mi##
## Moran I test under randomisation
##
## data: okanagan_sf$n_cyclist_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 14.984, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.4004592704 -0.0021505376 0.0007220019pd_cc_mi <- moran.test(okanagan_sf$n_pedestrian_casualty_claims, listw,zero.policy = TRUE)
pd_cc_mi##
## Moran I test under randomisation
##
## data: okanagan_sf$n_pedestrian_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 8.8431, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.2316536313 -0.0021505376 0.0006990321BYM2 Models
nb2INLA("okanagan.adj", okanagan_nb)
g <- inla.read.graph(filename = "okanagan.adj")
# Model Set 1: Total Casualty Claim Crashes
okanagan_models_1 <- cma_models(okanagan_sp@data, "Okanagan", "n_casualty_claims", "vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"okanagan_all <- clean_fixed_effects(okanagan_models_1$fixed_effects)
map(okanagan_models_1$models, ~ summary(.x))## $Unadjusted
## Time used:
## Pre = 14.1, Running = 0.829, Post = 0.166, Total = 15.1
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 2.164 0.047 2.072 2.164 2.255 2.164 0
## vandix_z 0.344 0.092 0.163 0.344 0.525 0.344 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.237 0.037 0.172 0.234 0.317 0.230
## Phi for ui 0.843 0.059 0.704 0.851 0.934 0.868
##
## Deviance Information Criterion (DIC) ...............: 2782.89
## Deviance Information Criterion (DIC, saturated) ....: 931.94
## Effective number of parameters .....................: 422.95
##
## Watanabe-Akaike information criterion (WAIC) ...: 2786.87
## Effective number of parameters .................: 293.18
##
## Marginal log-Likelihood: -1785.42
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.3, Running = 0.919, Post = 0.134, Total = 15.4
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 1.579 0.053 1.476 1.579 1.682 1.579 0
## vandix_z 0.464 0.074 0.318 0.463 0.610 0.463 0
## ln_roads_km_c 1.145 0.073 1.002 1.145 1.289 1.145 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.241 0.022 0.20 0.241 0.286 0.241
## Phi for ui 0.985 0.013 0.95 0.989 0.998 0.996
##
## Deviance Information Criterion (DIC) ...............: 2705.48
## Deviance Information Criterion (DIC, saturated) ....: 854.52
## Effective number of parameters .....................: 386.80
##
## Watanabe-Akaike information criterion (WAIC) ...: 2655.01
## Effective number of parameters .................: 237.31
##
## Marginal log-Likelihood: -1689.24
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.2, Running = 0.983, Post = 0.153, Total = 15.3
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode
## (Intercept) 1.962 0.084 1.796 1.962 2.126 1.962
## vandix_z 0.355 0.065 0.228 0.355 0.484 0.355
## ln_roads_km_c 0.949 0.076 0.801 0.949 1.098 0.949
## roads_prop_highway_arterial_z 0.605 0.054 0.499 0.605 0.711 0.605
## ale_index_z 0.510 0.263 -0.006 0.510 1.026 0.510
## canbics_index_z 0.016 0.215 -0.406 0.016 0.438 0.016
## population_100_c 0.031 0.009 0.014 0.031 0.048 0.031
## kld
## (Intercept) 0
## vandix_z 0
## ln_roads_km_c 0
## roads_prop_highway_arterial_z 0
## ale_index_z 0
## canbics_index_z 0
## population_100_c 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.339 0.036 0.271 0.339 0.412 0.340
## Phi for ui 0.979 0.019 0.929 0.985 0.998 0.995
##
## Deviance Information Criterion (DIC) ...............: 2695.73
## Deviance Information Criterion (DIC, saturated) ....: 844.77
## Effective number of parameters .....................: 371.48
##
## Watanabe-Akaike information criterion (WAIC) ...: 2644.06
## Effective number of parameters .................: 228.01
##
## Marginal log-Likelihood: -1650.75
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 2: Total Casualty Cyclist Claim Crashes
okanagan_models_2 <- cma_models(okanagan_sp@data,"Okanagan","n_cyclist_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"okanagan_cyclist <- clean_fixed_effects(okanagan_models_2$fixed_effects)
map(okanagan_models_2$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 14.1, Running = 0.748, Post = 0.142, Total = 15
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.113 0.110 -1.337 -1.111 -0.903 -1.112 0
## vandix_z 0.642 0.109 0.429 0.642 0.858 0.642 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.426 0.093 0.269 0.417 0.634 0.402
## Phi for ui 0.788 0.095 0.568 0.802 0.933 0.834
##
## Deviance Information Criterion (DIC) ...............: 995.71
## Deviance Information Criterion (DIC, saturated) ....: 526.98
## Effective number of parameters .....................: 160.92
##
## Watanabe-Akaike information criterion (WAIC) ...: 993.54
## Effective number of parameters .................: 117.40
##
## Marginal log-Likelihood: -488.26
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.1, Running = 0.832, Post = 0.155, Total = 15.1
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.507 0.140 -1.787 -1.505 -1.239 -1.505 0
## vandix_z 0.705 0.105 0.500 0.704 0.913 0.704 0
## ln_roads_km_c 0.721 0.121 0.484 0.721 0.961 0.721 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.339 0.062 0.23 0.334 0.474 0.327
## Phi for ui 0.951 0.041 0.84 0.963 0.995 0.985
##
## Deviance Information Criterion (DIC) ...............: 964.03
## Deviance Information Criterion (DIC, saturated) ....: 495.30
## Effective number of parameters .....................: 140.52
##
## Watanabe-Akaike information criterion (WAIC) ...: 957.29
## Effective number of parameters .................: 101.14
##
## Marginal log-Likelihood: -477.10
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.1, Running = 0.86, Post = 0.162, Total = 15.1
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -1.014 0.165 -1.346 -1.012 -0.697
## vandix_z 0.566 0.101 0.369 0.566 0.765
## ln_roads_km_c 0.616 0.137 0.348 0.616 0.887
## roads_prop_highway_arterial_z 0.395 0.089 0.221 0.394 0.569
## ale_index_z 0.703 0.379 -0.040 0.703 1.448
## canbics_index_z 0.354 0.292 -0.220 0.355 0.925
## population_100_c 0.036 0.015 0.008 0.036 0.065
## mode kld
## (Intercept) -1.012 0
## vandix_z 0.566 0
## ln_roads_km_c 0.616 0
## roads_prop_highway_arterial_z 0.394 0
## ale_index_z 0.703 0
## canbics_index_z 0.355 0
## population_100_c 0.036 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.428 0.082 0.288 0.421 0.608 0.409
## Phi for ui 0.961 0.038 0.858 0.973 0.997 0.991
##
## Deviance Information Criterion (DIC) ...............: 951.76
## Deviance Information Criterion (DIC, saturated) ....: 483.03
## Effective number of parameters .....................: 124.22
##
## Watanabe-Akaike information criterion (WAIC) ...: 945.65
## Effective number of parameters .................: 91.39
##
## Marginal log-Likelihood: -480.27
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 3: Total Casualty Cyclist Claim Crashes
okanagan_models_3 <- cma_models(okanagan_sp@data,"Okanagan","n_pedestrian_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"okanagan_pedestrian <- clean_fixed_effects(okanagan_models_3$fixed_effects)
map(okanagan_models_3$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 14.2, Running = 0.777, Post = 0.143, Total = 15.1
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.042 0.117 -1.284 -1.038 -0.822 -1.039 0
## vandix_z 0.668 0.110 0.452 0.668 0.884 0.668 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.525 0.088 0.374 0.517 0.720 0.501
## Phi for ui 0.328 0.123 0.121 0.317 0.593 0.291
##
## Deviance Information Criterion (DIC) ...............: 1114.17
## Deviance Information Criterion (DIC, saturated) ....: 620.03
## Effective number of parameters .....................: 220.56
##
## Watanabe-Akaike information criterion (WAIC) ...: 1132.39
## Effective number of parameters .................: 167.12
##
## Marginal log-Likelihood: -554.86
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.2, Running = 0.809, Post = 0.157, Total = 15.2
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.439 0.138 -1.716 -1.437 -1.173 -1.437 0
## vandix_z 0.731 0.119 0.497 0.730 0.965 0.730 0
## ln_roads_km_c 0.805 0.135 0.541 0.805 1.073 0.805 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.372 0.084 0.232 0.364 0.561 0.349
## Phi for ui 0.751 0.114 0.489 0.768 0.925 0.808
##
## Deviance Information Criterion (DIC) ...............: 1070.51
## Deviance Information Criterion (DIC, saturated) ....: 576.37
## Effective number of parameters .....................: 185.13
##
## Watanabe-Akaike information criterion (WAIC) ...: 1071.08
## Effective number of parameters .................: 135.15
##
## Marginal log-Likelihood: -542.90
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.4, Running = 0.914, Post = 0.149, Total = 15.4
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -0.861 0.160 -1.182 -0.859 -0.554
## vandix_z 0.558 0.113 0.338 0.558 0.780
## ln_roads_km_c 0.632 0.150 0.337 0.632 0.927
## roads_prop_highway_arterial_z 0.575 0.088 0.404 0.575 0.748
## ale_index_z 1.294 0.391 0.531 1.293 2.065
## canbics_index_z -0.245 0.321 -0.877 -0.244 0.384
## population_100_c 0.049 0.016 0.019 0.049 0.081
## mode kld
## (Intercept) -0.859 0
## vandix_z 0.558 0
## ln_roads_km_c 0.632 0
## roads_prop_highway_arterial_z 0.575 0
## ale_index_z 1.293 0
## canbics_index_z -0.244 0
## population_100_c 0.049 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.508 0.122 0.310 0.494 0.786 0.468
## Phi for ui 0.707 0.132 0.408 0.725 0.911 0.772
##
## Deviance Information Criterion (DIC) ...............: 1040.50
## Deviance Information Criterion (DIC, saturated) ....: 546.36
## Effective number of parameters .....................: 164.38
##
## Watanabe-Akaike information criterion (WAIC) ...: 1035.97
## Effective number of parameters .................: 118.86
##
## Marginal log-Likelihood: -534.59
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')okanagan_results <- bind_rows(okanagan_all %>% mutate(Region = "Okanagan",Outcome = "All Injury Claims"),
okanagan_cyclist %>% mutate(Region = "Okanagan",Outcome = "Cyclist Injury Claims"),
okanagan_pedestrian %>% mutate(Region = "Okanagan",Outcome = "Pedestrian Injury Claims")
) %>%
filter(variable == "vandix_z") %>%
select(Region,Outcome,everything())Fraser Valley
Descriptive Statistics
region | variable | n | sum | mean | sd | min | max | missing | n_zero | p_zero |
|---|---|---|---|---|---|---|---|---|---|---|
Fraser Valley | n_claims | 472 | 68,269.0 | 144.6 | 279.9 | 0.0 | 2,538.0 | 0 | 4 | 0.8 |
n_casualty_claims | 15,730.0 | 33.3 | 65.2 | 0.0 | 556.0 | 0 | 17 | 3.6 | ||
n_cyclist_claims | 581.0 | 1.2 | 2.4 | 0.0 | 21.0 | 0 | 248 | 52.5 | ||
n_cyclist_casualty_claims | 360.0 | 0.8 | 1.6 | 0.0 | 12.0 | 0 | 302 | 64.0 | ||
n_pedestrian_claims | 894.0 | 1.9 | 3.7 | 0.0 | 32.0 | 0 | 205 | 43.4 | ||
n_pedestrian_casualty_claims | 725.0 | 1.5 | 3.0 | 0.0 | 25.0 | 0 | 225 | 47.7 | ||
population | 309,493.0 | 655.7 | 341.8 | 0.0 | 2,550.0 | 0 | 4 | 0.8 | ||
total_roads_km | 3,678.9 | 7.8 | 9.6 | 0.1 | 81.8 | 0 | 0 | 0.0 | ||
roads_prop_highway_arterial | 40.5 | 0.1 | 0.1 | 0.0 | 0.7 | 0 | 281 | 59.5 | ||
no_highschool_prevalance | 82.9 | 0.2 | 0.1 | 0.1 | 0.5 | 88 | 0 | 0.0 | ||
unemployment_rate | 2,706.8 | 7.0 | 5.1 | 0.0 | 30.0 | 88 | 45 | 9.5 | ||
hh_avg_income | 27,748,698.0 | 75,404.1 | 21,091.5 | 24,502.0 | 161,236.0 | 104 | 0 | 0.0 | ||
participation_rate | 24,427.3 | 63.6 | 11.6 | 12.7 | 82.2 | 88 | 0 | 0.0 | ||
university_degree_prevalance | 170.5 | 0.4 | 0.1 | 0.1 | 0.7 | 88 | 0 | 0.0 | ||
lone_parent_fam_prevalence | 63.8 | 0.2 | 0.1 | 0.0 | 1.0 | 87 | 3 | 0.6 | ||
home_owner_prevalence | 284.5 | 0.7 | 0.2 | 0.0 | 1.0 | 88 | 1 | 0.2 | ||
vandix | 112.4 | 0.3 | 0.6 | -0.9 | 2.0 | 104 | 0 | 0.0 | ||
ale_index | -89.0 | -0.2 | 1.7 | -2.1 | 6.7 | 101 | 1 | 0.2 | ||
canbics_index | 514.5 | 1.4 | 1.3 | 0.0 | 4.9 | 101 | 81 | 17.2 |
Maps
claims <- ggplot() +
geom_sf(data = fraser_valley_sf, aes(fill = n_casualty_claims,colour=n_casualty_claims)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
theme_void() +
ggtitle(
"Vancouver"
)
vandix <- ggplot() +
geom_sf(data = fraser_valley_sf, aes(fill = vandix_z_c,colour=vandix_z_c)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
scale_colour_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
theme_void()
total_roads <- ggplot() +
geom_sf(data = fraser_valley_sf, aes(fill = total_roads_km,colour=total_roads_km)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
theme_void()
cowplot::plot_grid(claims,vandix,total_roads,ncol=1)
Setting up Spatial Weights Matrix
#### Define Spatial Neighrbourhoods
fraser_valley_sp <- as(fraser_valley_sf, "Spatial")
fraser_valley_sp$ui <- 1:nrow(fraser_valley_sp@data)
coords <- coordinates(fraser_valley_sp)
fraser_valley_nb <- poly2nb(fraser_valley_sp, queen = TRUE)## Warning in poly2nb(fraser_valley_sp, queen = TRUE): neighbour object has 2 sub-graphs;
## if this sub-graph count seems unexpected, try increasing the snap argument.fraser_valley_nb## Neighbour list object:
## Number of regions: 472
## Number of nonzero links: 2822
## Percentage nonzero weights: 1.266698
## Average number of links: 5.978814
## 2 disjoint connected subgraphs#assign nearest neighbour for no links
fraser_valley_nb <- assign_nearest_neighbors(fraser_valley_nb,fraser_valley_sp)
plot(fraser_valley_sp, border = grey(0.5))
plot(fraser_valley_nb,
coords = coords,
add = TRUE, pch = 16, lwd = 2)
Spatial Dependency in each crash type
listw <- nb2listw(fraser_valley_nb,zero.policy = TRUE)
all_cc_mi <- moran.test(fraser_valley_sf$n_casualty_claims, listw,zero.policy = TRUE)
all_cc_mi##
## Moran I test under randomisation
##
## data: fraser_valley_sf$n_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 8.2045, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.2194714130 -0.0021231423 0.0007294782cyc_cc_mi <- moran.test(fraser_valley_sf$n_cyclist_casualty_claims, listw,zero.policy = TRUE)
cyc_cc_mi##
## Moran I test under randomisation
##
## data: fraser_valley_sf$n_cyclist_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 8.2775, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.2236414065 -0.0021231423 0.0007438924pd_cc_mi <- moran.test(fraser_valley_sf$n_pedestrian_casualty_claims, listw,zero.policy = TRUE)
pd_cc_mi##
## Moran I test under randomisation
##
## data: fraser_valley_sf$n_pedestrian_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 9.3134, p-value < 2.2e-16
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.2508545992 -0.0021231423 0.0007378148BYM2 Models
nb2INLA("fraser_valley.adj", fraser_valley_nb)
g <- inla.read.graph(filename = "fraser_valley.adj")
# Model Set 1: Total Casualty Claim Crashes
fraser_valley_models_1 <- cma_models(fraser_valley_sp@data, "Fraser Valley", "n_casualty_claims", "vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"fraser_valley_all <- clean_fixed_effects(fraser_valley_models_1$fixed_effects)
map(fraser_valley_models_1$models, ~ summary(.x))## $Unadjusted
## Time used:
## Pre = 14.4, Running = 0.907, Post = 0.133, Total = 15.4
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 2.479 0.040 2.399 2.479 2.557 2.479 0
## vandix_z 0.232 0.085 0.065 0.232 0.399 0.232 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.299 0.022 0.258 0.298 0.344 0.297
## Phi for ui 0.767 0.354 0.001 0.985 1.000 1.000
##
## Deviance Information Criterion (DIC) ...............: 2987.34
## Deviance Information Criterion (DIC, saturated) ....: 924.07
## Effective number of parameters .....................: 430.45
##
## Watanabe-Akaike information criterion (WAIC) ...: 2950.25
## Effective number of parameters .................: 274.19
##
## Marginal log-Likelihood: -1818.76
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.4, Running = 0.84, Post = 0.14, Total = 15.4
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 2.265 0.043 2.180 2.265 2.350 2.265 0
## vandix_z 0.285 0.066 0.156 0.285 0.415 0.285 0
## ln_roads_km_c 1.104 0.067 0.972 1.104 1.236 1.104 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.712 0.069 0.585 0.709 0.856 0.703
## Phi for ui 0.802 0.055 0.680 0.807 0.893 0.818
##
## Deviance Information Criterion (DIC) ...............: 2985.16
## Deviance Information Criterion (DIC, saturated) ....: 921.90
## Effective number of parameters .....................: 417.80
##
## Watanabe-Akaike information criterion (WAIC) ...: 2950.49
## Effective number of parameters .................: 269.04
##
## Marginal log-Likelihood: -1674.82
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.5, Running = 0.949, Post = 0.184, Total = 15.6
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode
## (Intercept) 2.719 0.060 2.601 2.719 2.837 2.719
## vandix_z 0.140 0.060 0.022 0.139 0.257 0.139
## ln_roads_km_c 0.862 0.072 0.721 0.862 1.004 0.862
## roads_prop_highway_arterial_z 0.580 0.059 0.464 0.580 0.695 0.580
## ale_index_z 0.245 0.127 -0.004 0.245 0.494 0.245
## canbics_index_z 0.183 0.142 -0.097 0.183 0.458 0.183
## population_100_c 0.027 0.012 0.003 0.027 0.051 0.027
## kld
## (Intercept) 0
## vandix_z 0
## ln_roads_km_c 0
## roads_prop_highway_arterial_z 0
## ale_index_z 0
## canbics_index_z 0
## population_100_c 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.953 0.104 0.760 0.949 1.170 0.944
## Phi for ui 0.805 0.066 0.658 0.812 0.913 0.825
##
## Deviance Information Criterion (DIC) ...............: 2973.51
## Deviance Information Criterion (DIC, saturated) ....: 910.24
## Effective number of parameters .....................: 403.86
##
## Watanabe-Akaike information criterion (WAIC) ...: 2931.21
## Effective number of parameters .................: 256.46
##
## Marginal log-Likelihood: -1640.90
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 2: Total Casualty Cyclist Claim Crashes
fraser_valley_models_2 <- cma_models(fraser_valley_sp@data,"Fraser Valley","n_cyclist_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"fraser_valley_cyclist <- clean_fixed_effects(fraser_valley_models_2$fixed_effects)
map(fraser_valley_models_2$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 14.4, Running = 0.852, Post = 0.136, Total = 15.4
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.220 0.125 -1.478 -1.216 -0.984 -1.216 0
## vandix_z 0.338 0.113 0.113 0.338 0.557 0.338 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.679 0.113 0.485 0.669 0.930 0.648
## Phi for ui 0.315 0.130 0.098 0.303 0.593 0.275
##
## Deviance Information Criterion (DIC) ...............: 997.87
## Deviance Information Criterion (DIC, saturated) ....: 582.62
## Effective number of parameters .....................: 185.63
##
## Watanabe-Akaike information criterion (WAIC) ...: 1003.75
## Effective number of parameters .................: 138.73
##
## Marginal log-Likelihood: -315.14
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.4, Running = 1.03, Post = 0.189, Total = 15.6
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.360 0.124 -1.610 -1.358 -1.123 -1.358 0
## vandix_z 0.286 0.114 0.063 0.286 0.509 0.286 0
## ln_roads_km_c 0.753 0.127 0.505 0.752 1.004 0.752 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.623 0.114 0.427 0.614 0.874 0.596
## Phi for ui 0.752 0.120 0.479 0.769 0.934 0.813
##
## Deviance Information Criterion (DIC) ...............: 961.35
## Deviance Information Criterion (DIC, saturated) ....: 546.10
## Effective number of parameters .....................: 153.84
##
## Watanabe-Akaike information criterion (WAIC) ...: 964.81
## Effective number of parameters .................: 117.69
##
## Marginal log-Likelihood: -303.34
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.9, Running = 1.03, Post = 0.234, Total = 16.1
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -0.791 0.145 -1.082 -0.789 -0.513
## vandix_z 0.081 0.112 -0.139 0.081 0.299
## ln_roads_km_c 0.467 0.133 0.207 0.466 0.731
## roads_prop_highway_arterial_z 0.601 0.101 0.404 0.600 0.800
## ale_index_z 0.695 0.214 0.276 0.695 1.115
## canbics_index_z 0.091 0.299 -0.494 0.091 0.678
## population_100_c 0.046 0.023 0.001 0.046 0.091
## mode kld
## (Intercept) -0.789 0
## vandix_z 0.081 0
## ln_roads_km_c 0.466 0
## roads_prop_highway_arterial_z 0.600 0
## ale_index_z 0.695 0
## canbics_index_z 0.091 0
## population_100_c 0.046 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.897 0.186 0.588 0.879 1.316 0.842
## Phi for ui 0.611 0.160 0.278 0.623 0.882 0.659
##
## Deviance Information Criterion (DIC) ...............: 941.71
## Deviance Information Criterion (DIC, saturated) ....: 526.46
## Effective number of parameters .....................: 136.45
##
## Watanabe-Akaike information criterion (WAIC) ...: 942.40
## Effective number of parameters .................: 105.00
##
## Marginal log-Likelihood: -300.28
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 3: Total Casualty Cyclist Claim Crashes
fraser_valley_models_3 <- cma_models(fraser_valley_sp@data,"Fraser Valley","n_pedestrian_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"fraser_valley_pedestrian <- clean_fixed_effects(fraser_valley_models_3$fixed_effects)
map(fraser_valley_models_3$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 14.6, Running = 0.928, Post = 0.165, Total = 15.7
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -0.585 0.096 -0.779 -0.583 -0.402 -0.583 0
## vandix_z 0.503 0.095 0.318 0.503 0.689 0.503 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.740 0.104 0.558 0.732 0.965 0.716
## Phi for ui 0.372 0.115 0.166 0.366 0.610 0.356
##
## Deviance Information Criterion (DIC) ...............: 1347.02
## Deviance Information Criterion (DIC, saturated) ....: 695.41
## Effective number of parameters .....................: 236.66
##
## Watanabe-Akaike information criterion (WAIC) ...: 1353.70
## Effective number of parameters .................: 174.42
##
## Marginal log-Likelihood: -517.06
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.6, Running = 0.821, Post = 0.146, Total = 15.5
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -0.730 0.092 -0.915 -0.729 -0.552 -0.729 0
## vandix_z 0.447 0.091 0.269 0.447 0.627 0.447 0
## ln_roads_km_c 0.813 0.104 0.608 0.813 1.017 0.813 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.660 0.103 0.476 0.653 0.882 0.642
## Phi for ui 0.852 0.087 0.638 0.869 0.970 0.912
##
## Deviance Information Criterion (DIC) ...............: 1295.95
## Deviance Information Criterion (DIC, saturated) ....: 644.34
## Effective number of parameters .....................: 194.13
##
## Watanabe-Akaike information criterion (WAIC) ...: 1293.89
## Effective number of parameters .................: 143.01
##
## Marginal log-Likelihood: -495.02
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.5, Running = 0.87, Post = 0.169, Total = 15.5
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -0.210 0.114 -0.437 -0.208 0.010
## vandix_z 0.279 0.087 0.109 0.279 0.450
## ln_roads_km_c 0.421 0.114 0.199 0.421 0.646
## roads_prop_highway_arterial_z 0.619 0.082 0.459 0.619 0.781
## ale_index_z 0.488 0.178 0.140 0.488 0.837
## canbics_index_z 0.025 0.239 -0.445 0.026 0.494
## population_100_c 0.059 0.019 0.023 0.059 0.097
## mode kld
## (Intercept) -0.208 0
## vandix_z 0.279 0
## ln_roads_km_c 0.421 0
## roads_prop_highway_arterial_z 0.619 0
## ale_index_z 0.488 0
## canbics_index_z 0.026 0
## population_100_c 0.059 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1.044 0.180 0.732 1.030 1.440 1.002
## Phi for ui 0.675 0.131 0.392 0.688 0.892 0.719
##
## Deviance Information Criterion (DIC) ...............: 1271.74
## Deviance Information Criterion (DIC, saturated) ....: 620.13
## Effective number of parameters .....................: 173.99
##
## Watanabe-Akaike information criterion (WAIC) ...: 1264.84
## Effective number of parameters .................: 127.02
##
## Marginal log-Likelihood: -482.97
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')fraser_valley_results <- bind_rows(fraser_valley_all %>% mutate(Region = "Fraser Valley",Outcome = "All Injury Claims"),
fraser_valley_cyclist %>% mutate(Region = "Fraser Valley",Outcome = "Cyclist Injury Claims"),
fraser_valley_pedestrian %>% mutate(Region = "Fraser Valley",Outcome = "Pedestrian Injury Claims")
) %>%
filter(variable == "vandix_z") %>%
select(Region,Outcome,everything())North Central
Descriptive Statistics
region | variable | n | sum | mean | sd | min | max | missing | n_zero | p_zero |
|---|---|---|---|---|---|---|---|---|---|---|
North Central | n_claims | 359 | 31,422.0 | 87.5 | 180.5 | 0.0 | 1,906.0 | 0 | 22 | 6.1 |
n_casualty_claims | 4,490.0 | 12.5 | 25.2 | 0.0 | 228.0 | 0 | 57 | 15.9 | ||
n_cyclist_claims | 107.0 | 0.3 | 1.0 | 0.0 | 11.0 | 0 | 296 | 82.5 | ||
n_cyclist_casualty_claims | 54.0 | 0.2 | 0.5 | 0.0 | 5.0 | 0 | 323 | 90.0 | ||
n_pedestrian_claims | 241.0 | 0.7 | 2.0 | 0.0 | 27.0 | 0 | 249 | 69.4 | ||
n_pedestrian_casualty_claims | 203.0 | 0.6 | 1.6 | 0.0 | 21.0 | 0 | 262 | 73.0 | ||
population | 182,818.0 | 513.5 | 284.6 | 0.0 | 2,479.0 | 3 | 11 | 3.1 | ||
total_roads_km | 9,989.6 | 27.8 | 53.9 | 0.0 | 398.4 | 0 | 9 | 2.5 | ||
roads_prop_highway_arterial | 45.4 | 0.1 | 0.1 | 0.0 | 0.7 | 9 | 132 | 36.8 | ||
no_highschool_prevalance | 74.8 | 0.2 | 0.1 | 0.1 | 0.8 | 30 | 0 | 0.0 | ||
unemployment_rate | 3,459.5 | 10.5 | 6.8 | 0.0 | 60.0 | 30 | 14 | 3.9 | ||
hh_avg_income | 24,402,848.0 | 77,716.1 | 21,870.6 | 34,760.0 | 201,050.0 | 45 | 0 | 0.0 | ||
participation_rate | 22,327.9 | 67.9 | 9.3 | 35.7 | 89.8 | 30 | 0 | 0.0 | ||
university_degree_prevalance | 149.0 | 0.5 | 0.1 | 0.0 | 0.7 | 30 | 1 | 0.3 | ||
lone_parent_fam_prevalence | 60.1 | 0.2 | 0.1 | 0.0 | 0.6 | 29 | 0 | 0.0 | ||
home_owner_prevalence | 233.8 | 0.7 | 0.2 | 0.0 | 1.0 | 30 | 4 | 1.1 | ||
vandix | 156.6 | 0.5 | 0.7 | -1.0 | 3.4 | 45 | 0 | 0.0 | ||
ale_index | -387.8 | -1.3 | 0.6 | -2.1 | 0.0 | 65 | 0 | 0.0 | ||
canbics_index | 183.2 | 0.6 | 1.0 | 0.0 | 4.9 | 44 | 179 | 49.9 |
Maps
claims <- ggplot() +
geom_sf(data = north_central_sf, aes(fill = n_casualty_claims,colour=n_casualty_claims)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
theme_void() +
ggtitle(
"Vancouver"
)
vandix <- ggplot() +
geom_sf(data = north_central_sf, aes(fill = vandix_z_c,colour=vandix_z_c)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
scale_colour_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
theme_void()
total_roads <- ggplot() +
geom_sf(data = north_central_sf, aes(fill = total_roads_km,colour=total_roads_km)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
theme_void()
cowplot::plot_grid(claims,vandix,total_roads,ncol=1)
Setting up Spatial Weights Matrix
#### Define Spatial Neighrbourhoods
north_central_sp <- as(north_central_sf, "Spatial")
north_central_sp$ui <- 1:nrow(north_central_sp@data)
coords <- coordinates(north_central_sp)
north_central_nb <- poly2nb(north_central_sp, queen = TRUE)
north_central_nb## Neighbour list object:
## Number of regions: 359
## Number of nonzero links: 1958
## Percentage nonzero weights: 1.519231
## Average number of links: 5.454039#assign nearest neighbour for no links
north_central_nb <- assign_nearest_neighbors(north_central_nb,north_central_sp)
plot(north_central_sp, border = grey(0.5))
plot(north_central_nb,
coords = coords,
add = TRUE, pch = 16, lwd = 2)
Spatial Dependency in each crash type
listw <- nb2listw(north_central_nb,zero.policy = TRUE)
all_cc_mi <- moran.test(north_central_sf$n_casualty_claims, listw,zero.policy = TRUE)
all_cc_mi##
## Moran I test under randomisation
##
## data: north_central_sf$n_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 7.6708, p-value = 8.548e-15
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.241279277 -0.002793296 0.001012416cyc_cc_mi <- moran.test(north_central_sf$n_cyclist_casualty_claims, listw,zero.policy = TRUE)
cyc_cc_mi##
## Moran I test under randomisation
##
## data: north_central_sf$n_cyclist_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 3.7022, p-value = 0.0001069
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.116363461 -0.002793296 0.001035898pd_cc_mi <- moran.test(north_central_sf$n_pedestrian_casualty_claims, listw,zero.policy = TRUE)
pd_cc_mi##
## Moran I test under randomisation
##
## data: north_central_sf$n_pedestrian_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 6.5559, p-value = 2.766e-11
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.1938375579 -0.0027932961 0.0008995886BYM2 Models
nb2INLA("north_central.adj", north_central_nb)
g <- inla.read.graph(filename = "north_central.adj")
# Model Set 1: Total Casualty Claim Crashes
north_central_models_1 <- cma_models(north_central_sp@data, "North Central", "n_casualty_claims", "vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"north_central_all <- clean_fixed_effects(north_central_models_1$fixed_effects)
map(north_central_models_1$models, ~ summary(.x))## $Unadjusted
## Time used:
## Pre = 14.7, Running = 0.691, Post = 0.135, Total = 15.6
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 1.324 0.081 1.162 1.325 1.481 1.325 0
## vandix_z 0.191 0.076 0.042 0.191 0.340 0.191 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.296 0.044 0.218 0.293 0.391 0.287
## Phi for ui 0.788 0.062 0.650 0.794 0.890 0.807
##
## Deviance Information Criterion (DIC) ...............: 1843.35
## Deviance Information Criterion (DIC, saturated) ....: 690.46
## Effective number of parameters .....................: 301.98
##
## Watanabe-Akaike information criterion (WAIC) ...: 1868.79
## Effective number of parameters .................: 223.04
##
## Marginal log-Likelihood: -1102.66
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.2, Running = 0.685, Post = 0.13, Total = 15
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 0.531 0.136 0.261 0.533 0.795 0.533 0
## vandix_z 0.272 0.072 0.132 0.272 0.413 0.272 0
## ln_roads_km_c 0.698 0.090 0.521 0.698 0.876 0.698 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.229 0.031 0.172 0.228 0.294 0.227
## Phi for ui 0.950 0.031 0.871 0.957 0.990 0.972
##
## Deviance Information Criterion (DIC) ...............: 1801.00
## Deviance Information Criterion (DIC, saturated) ....: 648.12
## Effective number of parameters .....................: 283.39
##
## Watanabe-Akaike information criterion (WAIC) ...: 1797.13
## Effective number of parameters .................: 193.56
##
## Marginal log-Likelihood: -1079.96
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.1, Running = 0.755, Post = 0.156, Total = 15
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) 0.969 0.220 0.535 0.970 1.397
## vandix_z 0.201 0.067 0.071 0.201 0.332
## ln_roads_km_c 0.542 0.087 0.371 0.543 0.713
## roads_prop_highway_arterial_z 0.616 0.086 0.447 0.616 0.784
## ale_index_z -0.263 0.293 -0.837 -0.263 0.313
## canbics_index_z 0.384 0.311 -0.228 0.384 0.994
## population_100_c 0.034 0.021 -0.008 0.034 0.076
## mode kld
## (Intercept) 0.970 0
## vandix_z 0.201 0
## ln_roads_km_c 0.543 0
## roads_prop_highway_arterial_z 0.616 0
## ale_index_z -0.263 0
## canbics_index_z 0.384 0
## population_100_c 0.034 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.269 0.035 0.204 0.267 0.341 0.267
## Phi for ui 0.962 0.028 0.890 0.969 0.994 0.983
##
## Deviance Information Criterion (DIC) ...............: 1795.38
## Deviance Information Criterion (DIC, saturated) ....: 642.50
## Effective number of parameters .....................: 274.27
##
## Watanabe-Akaike information criterion (WAIC) ...: 1785.78
## Effective number of parameters .................: 185.29
##
## Marginal log-Likelihood: -1076.84
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 2: Total Casualty Cyclist Claim Crashes
north_central_models_2 <- cma_models(north_central_sp@data,"North Central","n_cyclist_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"north_central_cyclist <- clean_fixed_effects(north_central_models_2$fixed_effects)
map(north_central_models_2$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 14.1, Running = 0.734, Post = 0.119, Total = 15
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -3.130 0.345 -3.906 -3.11 -2.501 -3.117 0
## vandix_z 0.311 0.138 0.043 0.31 0.584 0.310 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.649 0.201 0.347 0.618 1.131 0.559
## Phi for ui 0.158 0.109 0.024 0.131 0.438 0.071
##
## Deviance Information Criterion (DIC) ...............: 277.05
## Deviance Information Criterion (DIC, saturated) ....: 195.66
## Effective number of parameters .....................: 53.64
##
## Watanabe-Akaike information criterion (WAIC) ...: 281.36
## Effective number of parameters .................: 43.50
##
## Marginal log-Likelihood: -75.31
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14, Running = 0.743, Post = 0.124, Total = 14.9
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -3.335 0.401 -4.194 -3.317 -2.595 -3.320 0
## vandix_z 0.350 0.145 0.069 0.349 0.637 0.349 0
## ln_roads_km_c 0.160 0.162 -0.154 0.159 0.483 0.159 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.637 0.198 0.340 0.607 1.113 0.549
## Phi for ui 0.199 0.129 0.032 0.169 0.522 0.097
##
## Deviance Information Criterion (DIC) ...............: 276.34
## Deviance Information Criterion (DIC, saturated) ....: 194.96
## Effective number of parameters .....................: 52.45
##
## Watanabe-Akaike information criterion (WAIC) ...: 280.80
## Effective number of parameters .................: 43.03
##
## Marginal log-Likelihood: -80.11
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.2, Running = 0.783, Post = 0.136, Total = 15.1
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -2.521 0.500 -3.542 -2.507 -1.578
## vandix_z 0.356 0.148 0.068 0.354 0.651
## ln_roads_km_c 0.211 0.181 -0.141 0.210 0.570
## roads_prop_highway_arterial_z 0.504 0.199 0.115 0.503 0.898
## ale_index_z -0.474 0.713 -1.872 -0.475 0.929
## canbics_index_z 2.028 0.574 0.911 2.024 3.168
## population_100_c -0.047 0.077 -0.197 -0.047 0.104
## mode kld
## (Intercept) -2.508 0
## vandix_z 0.354 0
## ln_roads_km_c 0.210 0
## roads_prop_highway_arterial_z 0.503 0
## ale_index_z -0.475 0
## canbics_index_z 2.024 0
## population_100_c -0.047 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.833 0.308 0.399 0.778 1.593 0.677
## Phi for ui 0.306 0.166 0.063 0.280 0.682 0.196
##
## Deviance Information Criterion (DIC) ...............: 272.87
## Deviance Information Criterion (DIC, saturated) ....: 191.48
## Effective number of parameters .....................: 43.18
##
## Watanabe-Akaike information criterion (WAIC) ...: 280.19
## Effective number of parameters .................: 39.40
##
## Marginal log-Likelihood: -89.42
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 3: Total Casualty Cyclist Claim Crashes
north_central_models_3 <- cma_models(north_central_sp@data,"North Central","n_pedestrian_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"north_central_pedestrian <- clean_fixed_effects(north_central_models_3$fixed_effects)
map(north_central_models_3$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 14.1, Running = 0.72, Post = 0.122, Total = 15
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.929 0.196 -2.344 -1.921 -1.566 -1.922 0
## vandix_z 0.382 0.105 0.177 0.382 0.589 0.382 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.548 0.124 0.346 0.534 0.833 0.506
## Phi for ui 0.359 0.152 0.106 0.346 0.681 0.308
##
## Deviance Information Criterion (DIC) ...............: 591.43
## Deviance Information Criterion (DIC, saturated) ....: 358.45
## Effective number of parameters .....................: 115.31
##
## Watanabe-Akaike information criterion (WAIC) ...: 593.68
## Effective number of parameters .................: 84.99
##
## Marginal log-Likelihood: -253.06
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.1, Running = 0.752, Post = 0.137, Total = 15
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -2.058 0.246 -2.566 -2.050 -1.598 -2.051 0
## vandix_z 0.392 0.107 0.182 0.392 0.604 0.392 0
## ln_roads_km_c 0.120 0.141 -0.147 0.116 0.409 0.116 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.512 0.132 0.302 0.495 0.818 0.463
## Phi for ui 0.451 0.180 0.131 0.445 0.801 0.425
##
## Deviance Information Criterion (DIC) ...............: 588.55
## Deviance Information Criterion (DIC, saturated) ....: 355.56
## Effective number of parameters .....................: 112.61
##
## Watanabe-Akaike information criterion (WAIC) ...: 590.64
## Effective number of parameters .................: 83.39
##
## Marginal log-Likelihood: -258.22
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.1, Running = 0.79, Post = 0.154, Total = 15.1
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -2.107 0.444 -3.013 -2.095 -1.271
## vandix_z 0.335 0.107 0.126 0.335 0.548
## ln_roads_km_c -0.063 0.162 -0.372 -0.066 0.265
## roads_prop_highway_arterial_z 0.650 0.146 0.364 0.650 0.939
## ale_index_z -1.238 0.631 -2.498 -1.230 -0.021
## canbics_index_z 0.642 0.503 -0.340 0.640 1.637
## population_100_c 0.022 0.044 -0.064 0.022 0.109
## mode kld
## (Intercept) -2.096 0
## vandix_z 0.335 0
## ln_roads_km_c -0.066 0
## roads_prop_highway_arterial_z 0.650 0
## ale_index_z -1.230 0
## canbics_index_z 0.640 0
## population_100_c 0.022 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.522 0.147 0.291 0.503 0.864 0.466
## Phi for ui 0.560 0.181 0.204 0.568 0.875 0.598
##
## Deviance Information Criterion (DIC) ...............: 575.77
## Deviance Information Criterion (DIC, saturated) ....: 342.79
## Effective number of parameters .....................: 103.09
##
## Watanabe-Akaike information criterion (WAIC) ...: 578.23
## Effective number of parameters .................: 77.80
##
## Marginal log-Likelihood: -266.22
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')north_central_results <- bind_rows(north_central_all %>% mutate(Region = "North Central",Outcome = "All Injury Claims"),
north_central_cyclist %>% mutate(Region = "North Central",Outcome = "Cyclist Injury Claims"),
north_central_pedestrian %>% mutate(Region = "North Central",Outcome = "Pedestrian Injury Claims")
) %>%
filter(variable == "vandix_z") %>%
select(Region,Outcome,everything())Kamloops-Salmon Arm
Descriptive Statistics
region | variable | n | sum | mean | sd | min | max | missing | n_zero | p_zero |
|---|---|---|---|---|---|---|---|---|---|---|
Kamloops-Salmon Arm | n_claims | 204 | 24,765.0 | 121.4 | 244.3 | 0.0 | 2,286.0 | 0 | 5 | 2.5 |
n_casualty_claims | 4,069.0 | 19.9 | 39.5 | 0.0 | 293.0 | 0 | 19 | 9.3 | ||
n_cyclist_claims | 143.0 | 0.7 | 1.7 | 0.0 | 14.0 | 0 | 138 | 67.6 | ||
n_cyclist_casualty_claims | 89.0 | 0.4 | 1.1 | 0.0 | 10.0 | 0 | 152 | 74.5 | ||
n_pedestrian_claims | 236.0 | 1.2 | 2.9 | 0.0 | 21.0 | 0 | 130 | 63.7 | ||
n_pedestrian_casualty_claims | 183.0 | 0.9 | 2.3 | 0.0 | 17.0 | 0 | 137 | 67.2 | ||
population | 133,847.0 | 656.1 | 348.0 | 0.0 | 2,310.0 | 0 | 2 | 1.0 | ||
total_roads_km | 2,744.5 | 13.5 | 25.3 | 0.5 | 217.0 | 0 | 0 | 0.0 | ||
roads_prop_highway_arterial | 49.9 | 0.2 | 0.2 | 0.0 | 0.7 | 0 | 32 | 15.7 | ||
no_highschool_prevalance | 29.4 | 0.2 | 0.1 | 0.0 | 0.4 | 31 | 0 | 0.0 | ||
unemployment_rate | 1,439.0 | 8.3 | 5.3 | 0.0 | 37.5 | 31 | 12 | 5.9 | ||
hh_avg_income | 12,008,631.0 | 71,479.9 | 20,375.9 | 29,162.0 | 137,275.0 | 36 | 0 | 0.0 | ||
participation_rate | 10,775.2 | 62.3 | 9.9 | 30.6 | 90.0 | 31 | 0 | 0.0 | ||
university_degree_prevalance | 88.4 | 0.5 | 0.1 | 0.3 | 0.7 | 31 | 0 | 0.0 | ||
lone_parent_fam_prevalence | 31.0 | 0.2 | 0.1 | 0.0 | 0.5 | 30 | 0 | 0.0 | ||
home_owner_prevalence | 128.1 | 0.7 | 0.2 | 0.1 | 1.0 | 31 | 0 | 0.0 | ||
vandix | 30.5 | 0.2 | 0.7 | -1.3 | 2.4 | 36 | 0 | 0.0 | ||
ale_index | -116.3 | -0.7 | 1.2 | -2.1 | 3.3 | 41 | 1 | 0.5 | ||
canbics_index | 227.8 | 1.4 | 1.5 | 0.0 | 5.6 | 41 | 40 | 19.6 |
Maps
claims <- ggplot() +
geom_sf(data = interior_sf, aes(fill = n_casualty_claims,colour=n_casualty_claims)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
theme_void() +
ggtitle(
"Vancouver"
)
vandix <- ggplot() +
geom_sf(data = interior_sf, aes(fill = vandix_z_c,colour=vandix_z_c)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
scale_colour_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
theme_void()
total_roads <- ggplot() +
geom_sf(data = interior_sf, aes(fill = total_roads_km,colour=total_roads_km)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
theme_void()
cowplot::plot_grid(claims,vandix,total_roads,ncol=1)
Setting up Spatial Weights Matrix
#### Define Spatial Neighrbourhoods
interior_sp <- as(interior_sf, "Spatial")
interior_sp$ui <- 1:nrow(interior_sp@data)
coords <- coordinates(interior_sp)
interior_nb <- poly2nb(interior_sp, queen = TRUE)## Warning in poly2nb(interior_sp, queen = TRUE): neighbour object has 3 sub-graphs;
## if this sub-graph count seems unexpected, try increasing the snap argument.interior_nb## Neighbour list object:
## Number of regions: 204
## Number of nonzero links: 1098
## Percentage nonzero weights: 2.638408
## Average number of links: 5.382353
## 3 disjoint connected subgraphs#assign nearest neighbour for no links
interior_nb <- assign_nearest_neighbors(interior_nb,interior_sp)
plot(interior_sp, border = grey(0.5))
plot(interior_nb,
coords = coords,
add = TRUE, pch = 16, lwd = 2)
Spatial Dependency in each crash type
listw <- nb2listw(interior_nb,zero.policy = TRUE)
all_cc_mi <- moran.test(interior_sf$n_casualty_claims, listw,zero.policy = TRUE)
all_cc_mi##
## Moran I test under randomisation
##
## data: interior_sf$n_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 6.3114, p-value = 1.383e-10
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.261379586 -0.004926108 0.001780360cyc_cc_mi <- moran.test(interior_sf$n_cyclist_casualty_claims, listw,zero.policy = TRUE)
cyc_cc_mi##
## Moran I test under randomisation
##
## data: interior_sf$n_cyclist_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 2.4815, p-value = 0.006542
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.096406801 -0.004926108 0.001667545pd_cc_mi <- moran.test(interior_sf$n_pedestrian_casualty_claims, listw,zero.policy = TRUE)
pd_cc_mi##
## Moran I test under randomisation
##
## data: interior_sf$n_pedestrian_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 4.6131, p-value = 1.983e-06
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.190275369 -0.004926108 0.001790493BYM2 Models
nb2INLA("interior.adj", interior_nb)
g <- inla.read.graph(filename = "interior.adj")
# Model Set 1: Total Casualty Claim Crashes
interior_models_1 <- cma_models(interior_sp@data, "Kamloops-Salmon Arm", "n_casualty_claims", "vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"interior_all <- clean_fixed_effects(interior_models_1$fixed_effects)
map(interior_models_1$models, ~ summary(.x))## $Unadjusted
## Time used:
## Pre = 14.6, Running = 0.538, Post = 0.114, Total = 15.3
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 1.953 0.082 1.790 1.954 2.112 1.954 0
## vandix_z 0.250 0.115 0.024 0.251 0.476 0.251 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.647 0.083 0.498 0.642 0.825 0.632
## Phi for ui 0.421 0.107 0.223 0.417 0.639 0.409
##
## Deviance Information Criterion (DIC) ...............: 1176.14
## Deviance Information Criterion (DIC, saturated) ....: 411.35
## Effective number of parameters .....................: 182.29
##
## Watanabe-Akaike information criterion (WAIC) ...: 1181.66
## Effective number of parameters .................: 129.25
##
## Marginal log-Likelihood: -647.01
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.6, Running = 0.519, Post = 0.127, Total = 15.3
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 1.431 0.099 1.234 1.431 1.624 1.431 0
## vandix_z 0.258 0.104 0.055 0.258 0.461 0.258 0
## ln_roads_km_c 0.850 0.106 0.644 0.849 1.058 0.849 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.776 0.104 0.590 0.77 0.999 0.759
## Phi for ui 0.625 0.097 0.425 0.63 0.802 0.637
##
## Deviance Information Criterion (DIC) ...............: 1158.35
## Deviance Information Criterion (DIC, saturated) ....: 393.55
## Effective number of parameters .....................: 172.37
##
## Watanabe-Akaike information criterion (WAIC) ...: 1151.76
## Effective number of parameters .................: 115.90
##
## Marginal log-Likelihood: -622.71
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.5, Running = 0.553, Post = 0.145, Total = 15.2
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) 1.132 0.119 0.898 1.133 1.366
## vandix_z 0.208 0.085 0.042 0.208 0.375
## ln_roads_km_c 0.491 0.096 0.302 0.491 0.681
## roads_prop_highway_arterial_z 0.534 0.056 0.424 0.534 0.645
## ale_index_z 0.091 0.257 -0.413 0.091 0.594
## canbics_index_z -0.712 0.249 -1.202 -0.712 -0.224
## population_100_c 0.099 0.019 0.062 0.099 0.136
## mode kld
## (Intercept) 1.133 0
## vandix_z 0.208 0
## ln_roads_km_c 0.491 0
## roads_prop_highway_arterial_z 0.534 0
## ale_index_z 0.091 0
## canbics_index_z -0.712 0
## population_100_c 0.099 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1.176 0.170 0.871 1.166 1.540 1.149
## Phi for ui 0.819 0.085 0.620 0.832 0.948 0.862
##
## Deviance Information Criterion (DIC) ...............: 1138.16
## Deviance Information Criterion (DIC, saturated) ....: 373.37
## Effective number of parameters .....................: 154.77
##
## Watanabe-Akaike information criterion (WAIC) ...: 1123.45
## Effective number of parameters .................: 100.52
##
## Marginal log-Likelihood: -598.04
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 2: Total Casualty Cyclist Claim Crashes
interior_models_2 <- cma_models(interior_sp@data,"Kamloops-Salmon Arm","n_cyclist_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"interior_cyclist <- clean_fixed_effects(interior_models_2$fixed_effects)
map(interior_models_2$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 14.5, Running = 0.524, Post = 0.131, Total = 15.2
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.615 0.208 -2.054 -1.608 -1.228 -1.609 0
## vandix_z 0.372 0.154 0.064 0.373 0.668 0.374 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1.061 0.334 0.563 1.009 1.863 0.911
## Phi for ui 0.209 0.131 0.034 0.182 0.528 0.106
##
## Deviance Information Criterion (DIC) ...............: 323.23
## Deviance Information Criterion (DIC, saturated) ....: 202.50
## Effective number of parameters .....................: 50.21
##
## Watanabe-Akaike information criterion (WAIC) ...: 325.30
## Effective number of parameters .................: 40.32
##
## Marginal log-Likelihood: -58.43
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.5, Running = 0.53, Post = 0.147, Total = 15.1
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -2.040 0.278 -2.618 -2.030 -1.523 -2.031 0
## vandix_z 0.394 0.163 0.069 0.395 0.709 0.395 0
## ln_roads_km_c 0.586 0.201 0.200 0.583 0.990 0.583 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1.066 0.349 0.551 1.01 1.909 0.906
## Phi for ui 0.519 0.193 0.160 0.52 0.871 0.514
##
## Deviance Information Criterion (DIC) ...............: 316.59
## Deviance Information Criterion (DIC, saturated) ....: 195.85
## Effective number of parameters .....................: 45.63
##
## Watanabe-Akaike information criterion (WAIC) ...: 318.57
## Effective number of parameters .................: 37.17
##
## Marginal log-Likelihood: -59.22
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.5, Running = 0.568, Post = 0.122, Total = 15.2
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -1.919 0.324 -2.590 -1.907 -1.318
## vandix_z 0.391 0.177 0.034 0.395 0.726
## ln_roads_km_c 0.195 0.216 -0.224 0.193 0.624
## roads_prop_highway_arterial_z 0.395 0.138 0.125 0.394 0.668
## ale_index_z 0.630 0.488 -0.340 0.635 1.576
## canbics_index_z -0.284 0.489 -1.239 -0.286 0.683
## population_100_c 0.141 0.038 0.067 0.141 0.217
## mode kld
## (Intercept) -1.907 0
## vandix_z 0.395 0
## ln_roads_km_c 0.193 0
## roads_prop_highway_arterial_z 0.394 0
## ale_index_z 0.635 0
## canbics_index_z -0.286 0
## population_100_c 0.141 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1.847 0.897 0.722 1.648 4.15 1.324
## Phi for ui 0.619 0.218 0.170 0.646 0.95 0.803
##
## Deviance Information Criterion (DIC) ...............: 308.25
## Deviance Information Criterion (DIC, saturated) ....: 187.51
## Effective number of parameters .....................: 34.28
##
## Watanabe-Akaike information criterion (WAIC) ...: 311.20
## Effective number of parameters .................: 30.41
##
## Marginal log-Likelihood: -66.20
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 3: Total Casualty Cyclist Claim Crashes
interior_models_3 <- cma_models(interior_sp@data,"Kamloops-Salmon Arm","n_pedestrian_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"interior_pedestrian <- clean_fixed_effects(interior_models_3$fixed_effects)
map(interior_models_3$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 14.5, Running = 0.532, Post = 0.116, Total = 15.1
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.463 0.201 -1.886 -1.455 -1.091 -1.456 0
## vandix_z 0.654 0.159 0.343 0.653 0.970 0.653 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.626 0.133 0.404 0.612 0.925 0.586
## Phi for ui 0.267 0.124 0.075 0.250 0.549 0.207
##
## Deviance Information Criterion (DIC) ...............: 415.59
## Deviance Information Criterion (DIC, saturated) ....: 244.14
## Effective number of parameters .....................: 81.09
##
## Watanabe-Akaike information criterion (WAIC) ...: 423.72
## Effective number of parameters .................: 63.39
##
## Marginal log-Likelihood: -121.84
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.6, Running = 0.53, Post = 0.149, Total = 15.2
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.851 0.243 -2.350 -1.843 -1.394 -1.844 0
## vandix_z 0.693 0.162 0.377 0.692 1.012 0.692 0
## ln_roads_km_c 0.620 0.181 0.269 0.619 0.980 0.619 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.676 0.152 0.426 0.66 1.020 0.628
## Phi for ui 0.474 0.155 0.190 0.47 0.778 0.454
##
## Deviance Information Criterion (DIC) ...............: 408.87
## Deviance Information Criterion (DIC, saturated) ....: 237.41
## Effective number of parameters .....................: 73.21
##
## Watanabe-Akaike information criterion (WAIC) ...: 414.87
## Effective number of parameters .................: 57.58
##
## Marginal log-Likelihood: -121.42
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.5, Running = 0.564, Post = 0.123, Total = 15.2
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -1.579 0.290 -2.173 -1.570 -1.035
## vandix_z 0.665 0.153 0.362 0.665 0.966
## ln_roads_km_c 0.248 0.189 -0.120 0.247 0.622
## roads_prop_highway_arterial_z 0.389 0.123 0.148 0.388 0.633
## ale_index_z 1.566 0.450 0.675 1.568 2.446
## canbics_index_z -0.983 0.452 -1.881 -0.980 -0.103
## population_100_c 0.111 0.037 0.039 0.111 0.184
## mode kld
## (Intercept) -1.570 0
## vandix_z 0.665 0
## ln_roads_km_c 0.247 0
## roads_prop_highway_arterial_z 0.388 0
## ale_index_z 1.568 0
## canbics_index_z -0.980 0
## population_100_c 0.111 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1.083 0.313 0.602 1.038 1.823 0.953
## Phi for ui 0.373 0.197 0.064 0.352 0.783 0.253
##
## Deviance Information Criterion (DIC) ...............: 400.81
## Deviance Information Criterion (DIC, saturated) ....: 229.36
## Effective number of parameters .....................: 61.39
##
## Watanabe-Akaike information criterion (WAIC) ...: 404.56
## Effective number of parameters .................: 49.31
##
## Marginal log-Likelihood: -125.78
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')interior_results <- bind_rows(interior_all %>% mutate(Region = "Kamloops-Salmon Arm",Outcome = "All Injury Claims"),
interior_cyclist %>% mutate(Region = "Kamloops-Salmon Arm",Outcome = "Cyclist Injury Claims"),
interior_pedestrian %>% mutate(Region = "Kamloops-Salmon Arm",Outcome = "Pedestrian Injury Claims")
) %>%
filter(variable == "vandix_z") %>%
select(Region,Outcome,everything())Southeast
Descriptive Statistics
region | variable | n | sum | mean | sd | min | max | missing | n_zero | p_zero |
|---|---|---|---|---|---|---|---|---|---|---|
Southeast | n_claims | 110 | 9,978.0 | 90.7 | 170.1 | 0.0 | 1,235.0 | 0 | 5 | 4.5 |
n_casualty_claims | 1,037.0 | 9.4 | 17.8 | 0.0 | 121.0 | 0 | 14 | 12.7 | ||
n_cyclist_claims | 40.0 | 0.4 | 0.8 | 0.0 | 4.0 | 0 | 87 | 79.1 | ||
n_cyclist_casualty_claims | 24.0 | 0.2 | 0.6 | 0.0 | 4.0 | 0 | 94 | 85.5 | ||
n_pedestrian_claims | 62.0 | 0.6 | 1.3 | 0.0 | 8.0 | 0 | 79 | 71.8 | ||
n_pedestrian_casualty_claims | 51.0 | 0.5 | 1.2 | 0.0 | 8.0 | 0 | 85 | 77.3 | ||
population | 60,427.0 | 549.3 | 214.3 | 0.0 | 1,482.0 | 0 | 4 | 3.6 | ||
total_roads_km | 1,598.6 | 14.5 | 23.7 | 0.0 | 202.4 | 0 | 2 | 1.8 | ||
roads_prop_highway_arterial | 13.9 | 0.1 | 0.2 | 0.0 | 0.5 | 2 | 49 | 44.5 | ||
no_highschool_prevalance | 16.2 | 0.2 | 0.1 | 0.0 | 0.3 | 6 | 0 | 0.0 | ||
unemployment_rate | 838.5 | 8.1 | 4.4 | 0.0 | 17.5 | 6 | 9 | 8.2 | ||
hh_avg_income | 7,149,722.0 | 68,747.3 | 16,755.0 | 34,159.0 | 122,980.0 | 6 | 0 | 0.0 | ||
participation_rate | 6,473.1 | 62.2 | 8.5 | 40.9 | 80.0 | 6 | 0 | 0.0 | ||
university_degree_prevalance | 57.1 | 0.5 | 0.1 | 0.3 | 0.7 | 6 | 0 | 0.0 | ||
lone_parent_fam_prevalence | 17.1 | 0.2 | 0.1 | 0.0 | 0.4 | 6 | 0 | 0.0 | ||
home_owner_prevalence | 77.2 | 0.7 | 0.2 | 0.0 | 1.0 | 6 | 1 | 0.9 | ||
vandix | 8.2 | 0.1 | 0.6 | -0.9 | 1.6 | 6 | 0 | 0.0 | ||
ale_index | -125.7 | -1.2 | 0.7 | -2.1 | 0.4 | 9 | 0 | 0.0 | ||
canbics_index | 32.4 | 0.3 | 0.5 | 0.0 | 2.0 | 9 | 64 | 58.2 |
Maps
claims <- ggplot() +
geom_sf(data = southeast_sf, aes(fill = n_casualty_claims,colour=n_casualty_claims)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
theme_void() +
ggtitle(
"Vancouver"
)
vandix <- ggplot() +
geom_sf(data = southeast_sf, aes(fill = vandix_z_c,colour=vandix_z_c)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
scale_colour_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
theme_void()
total_roads <- ggplot() +
geom_sf(data = southeast_sf, aes(fill = total_roads_km,colour=total_roads_km)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
theme_void()
cowplot::plot_grid(claims,vandix,total_roads,ncol=1)
Setting up Spatial Weights Matrix
#### Define Spatial Neighrbourhoods
southeast_sp <- as(southeast_sf, "Spatial")
southeast_sp$ui <- 1:nrow(southeast_sp@data)
coords <- coordinates(southeast_sp)
southeast_nb <- poly2nb(southeast_sp, queen = TRUE)## Warning in poly2nb(southeast_sp, queen = TRUE): neighbour object has 3 sub-graphs;
## if this sub-graph count seems unexpected, try increasing the snap argument.southeast_nb## Neighbour list object:
## Number of regions: 110
## Number of nonzero links: 552
## Percentage nonzero weights: 4.561983
## Average number of links: 5.018182
## 3 disjoint connected subgraphs#assign nearest neighbour for no links
southeast_nb <- assign_nearest_neighbors(southeast_nb,southeast_sp)
plot(southeast_sp, border = grey(0.5))
plot(southeast_nb,
coords = coords,
add = TRUE, pch = 16, lwd = 2)
Spatial Dependency in each crash type
listw <- nb2listw(southeast_nb,zero.policy = TRUE)
all_cc_mi <- moran.test(southeast_sf$n_casualty_claims, listw,zero.policy = TRUE)
all_cc_mi##
## Moran I test under randomisation
##
## data: southeast_sf$n_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 2.7069, p-value = 0.003396
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.146644926 -0.009174312 0.003313622cyc_cc_mi <- moran.test(southeast_sf$n_cyclist_casualty_claims, listw,zero.policy = TRUE)
cyc_cc_mi##
## Moran I test under randomisation
##
## data: southeast_sf$n_cyclist_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 0.65675, p-value = 0.2557
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.029151043 -0.009174312 0.003405459pd_cc_mi <- moran.test(southeast_sf$n_pedestrian_casualty_claims, listw,zero.policy = TRUE)
pd_cc_mi##
## Moran I test under randomisation
##
## data: southeast_sf$n_pedestrian_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 2.1446, p-value = 0.01599
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.114119689 -0.009174312 0.003305202BYM2 Models
nb2INLA("southeast.adj", southeast_nb)
g <- inla.read.graph(filename = "southeast.adj")
# Model Set 1: Total Casualty Claim Crashes
southeast_models_1 <- cma_models(southeast_sp@data, "Southeast", "n_casualty_claims", "vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"southeast_all <- clean_fixed_effects(southeast_models_1$fixed_effects)
map(southeast_models_1$models, ~ summary(.x))## $Unadjusted
## Time used:
## Pre = 14.4, Running = 0.42, Post = 0.14, Total = 14.9
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 1.270 0.111 1.046 1.271 1.484 1.271 0
## vandix_z 0.352 0.150 0.057 0.352 0.648 0.352 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.750 0.128 0.528 0.740 1.029 0.721
## Phi for ui 0.396 0.141 0.151 0.387 0.689 0.361
##
## Deviance Information Criterion (DIC) ...............: 543.48
## Deviance Information Criterion (DIC, saturated) ....: 207.17
## Effective number of parameters .....................: 89.95
##
## Watanabe-Akaike information criterion (WAIC) ...: 547.74
## Effective number of parameters .................: 64.68
##
## Marginal log-Likelihood: -266.28
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.4, Running = 0.434, Post = 0.15, Total = 14.9
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 0.651 0.166 0.321 0.652 0.973 0.652 0
## vandix_z 0.386 0.141 0.111 0.385 0.664 0.385 0
## ln_roads_km_c 0.779 0.157 0.474 0.778 1.092 0.778 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.827 0.147 0.571 0.815 1.150 0.795
## Phi for ui 0.691 0.145 0.377 0.707 0.923 0.754
##
## Deviance Information Criterion (DIC) ...............: 530.90
## Deviance Information Criterion (DIC, saturated) ....: 194.60
## Effective number of parameters .....................: 83.84
##
## Watanabe-Akaike information criterion (WAIC) ...: 528.31
## Effective number of parameters .................: 56.64
##
## Marginal log-Likelihood: -259.52
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.5, Running = 0.524, Post = 0.126, Total = 15.1
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) 1.564 0.397 0.776 1.566 2.338
## vandix_z 0.350 0.131 0.094 0.349 0.609
## ln_roads_km_c 0.490 0.179 0.142 0.489 0.844
## roads_prop_highway_arterial_z 0.627 0.150 0.331 0.627 0.921
## ale_index_z 1.634 0.820 0.032 1.631 3.253
## canbics_index_z -0.712 0.851 -2.390 -0.711 0.955
## population_100_c 0.152 0.057 0.040 0.151 0.265
## mode kld
## (Intercept) 1.566 0
## vandix_z 0.349 0
## ln_roads_km_c 0.489 0
## roads_prop_highway_arterial_z 0.627 0
## ale_index_z 1.631 0
## canbics_index_z -0.711 0
## population_100_c 0.151 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.994 0.183 0.677 0.979 1.397 0.952
## Phi for ui 0.771 0.134 0.457 0.795 0.961 0.870
##
## Deviance Information Criterion (DIC) ...............: 524.87
## Deviance Information Criterion (DIC, saturated) ....: 188.56
## Effective number of parameters .....................: 79.56
##
## Watanabe-Akaike information criterion (WAIC) ...: 518.83
## Effective number of parameters .................: 52.02
##
## Marginal log-Likelihood: -267.01
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 2: Total Casualty Cyclist Claim Crashes
southeast_models_2 <- cma_models(southeast_sp@data,"Southeast","n_cyclist_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"southeast_cyclist <- clean_fixed_effects(southeast_models_2$fixed_effects)
map(southeast_models_2$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 14.4, Running = 0.4, Post = 0.147, Total = 15
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.997 0.353 -2.783 -1.966 -1.387 -1.913 0
## vandix_z 0.597 0.230 0.156 0.593 1.061 0.593 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 11.145 30.615 0.759 4.049 67.116 1.404
## Phi for ui 0.277 0.217 0.018 0.218 0.796 0.048
##
## Deviance Information Criterion (DIC) ...............: 125.20
## Deviance Information Criterion (DIC, saturated) ....: 88.73
## Effective number of parameters .....................: 9.52
##
## Watanabe-Akaike information criterion (WAIC) ...: 133.43
## Effective number of parameters .................: 15.36
##
## Marginal log-Likelihood: 11.87
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.4, Running = 0.424, Post = 0.126, Total = 14.9
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -2.450 0.436 -3.362 -2.433 -1.643 -2.434 0
## vandix_z 0.755 0.249 0.271 0.753 1.249 0.753 0
## ln_roads_km_c 0.468 0.243 -0.008 0.467 0.948 0.467 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 60.416 409.765 0.709 7.817 408.21 1.435
## Phi for ui 0.357 0.246 0.027 0.308 0.88 0.078
##
## Deviance Information Criterion (DIC) ...............: 123.94
## Deviance Information Criterion (DIC, saturated) ....: 87.47
## Effective number of parameters .....................: 9.67
##
## Watanabe-Akaike information criterion (WAIC) ...: 128.56
## Effective number of parameters .................: 12.56
##
## Marginal log-Likelihood: 9.28
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.5, Running = 0.425, Post = 0.128, Total = 15.1
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -0.882 0.610 -2.094 -0.877 0.301
## vandix_z 0.828 0.252 0.337 0.828 1.325
## ln_roads_km_c 0.202 0.376 -0.535 0.202 0.939
## roads_prop_highway_arterial_z 0.331 0.334 -0.323 0.330 0.989
## ale_index_z 1.940 1.360 -0.715 1.936 4.622
## canbics_index_z 0.074 1.281 -2.441 0.075 2.586
## population_100_c 0.288 0.094 0.105 0.288 0.473
## mode kld
## (Intercept) -0.878 0
## vandix_z 0.828 0
## ln_roads_km_c 0.202 0
## roads_prop_highway_arterial_z 0.330 0
## ale_index_z 1.937 0
## canbics_index_z 0.075 0
## population_100_c 0.288 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 5986.857 1.50e+05 0.933 66.570 2.28e+04 0.970
## Phi for ui 0.353 2.63e-01 0.018 0.293 9.06e-01 0.041
##
## Deviance Information Criterion (DIC) ...............: 120.86
## Deviance Information Criterion (DIC, saturated) ....: 84.39
## Effective number of parameters .....................: 7.95
##
## Watanabe-Akaike information criterion (WAIC) ...: 124.65
## Effective number of parameters .................: 10.26
##
## Marginal log-Likelihood: -0.811
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 3: Total Casualty Cyclist Claim Crashes
southeast_models_3 <- cma_models(southeast_sp@data,"Southeast","n_pedestrian_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"southeast_pedestrian <- clean_fixed_effects(southeast_models_3$fixed_effects)
map(southeast_models_3$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 14.4, Running = 0.443, Post = 0.112, Total = 14.9
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.868 0.318 -2.559 -1.851 -1.289 -1.807 0
## vandix_z 0.545 0.240 0.075 0.543 1.021 0.543 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.728 0.242 0.368 0.689 1.310 0.618
## Phi for ui 0.165 0.143 0.011 0.122 0.546 0.031
##
## Deviance Information Criterion (DIC) ...............: 167.70
## Deviance Information Criterion (DIC, saturated) ....: 107.52
## Effective number of parameters .....................: 33.05
##
## Watanabe-Akaike information criterion (WAIC) ...: 171.01
## Effective number of parameters .................: 26.27
##
## Marginal log-Likelihood: -23.60
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14.5, Running = 0.408, Post = 0.157, Total = 15
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -2.101 0.378 -2.898 -2.085 -1.403 -2.087 0
## vandix_z 0.597 0.247 0.113 0.597 1.085 0.597 0
## ln_roads_km_c 0.295 0.257 -0.211 0.295 0.801 0.295 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.753 0.260 0.372 0.711 1.382 0.633
## Phi for ui 0.204 0.167 0.014 0.156 0.634 0.039
##
## Deviance Information Criterion (DIC) ...............: 168.93
## Deviance Information Criterion (DIC, saturated) ....: 108.74
## Effective number of parameters .....................: 32.75
##
## Watanabe-Akaike information criterion (WAIC) ...: 173.11
## Effective number of parameters .................: 26.76
##
## Marginal log-Likelihood: -27.78
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14.6, Running = 0.475, Post = 0.147, Total = 15.2
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -1.384 0.793 -2.998 -1.366 0.125
## vandix_z 0.606 0.253 0.111 0.605 1.106
## ln_roads_km_c -0.070 0.387 -0.840 -0.067 0.685
## roads_prop_highway_arterial_z 0.307 0.332 -0.349 0.308 0.958
## ale_index_z 0.615 1.423 -2.208 0.622 3.394
## canbics_index_z -0.167 1.444 -2.983 -0.177 2.704
## population_100_c 0.233 0.115 0.010 0.232 0.464
## mode kld
## (Intercept) -1.367 0
## vandix_z 0.606 0
## ln_roads_km_c -0.067 0
## roads_prop_highway_arterial_z 0.308 0
## ale_index_z 0.622 0
## canbics_index_z -0.177 0
## population_100_c 0.232 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.736 0.269 0.349 0.690 1.392 0.606
## Phi for ui 0.227 0.180 0.016 0.177 0.679 0.045
##
## Deviance Information Criterion (DIC) ...............: 170.16
## Deviance Information Criterion (DIC, saturated) ....: 109.97
## Effective number of parameters .....................: 33.93
##
## Watanabe-Akaike information criterion (WAIC) ...: 177.67
## Effective number of parameters .................: 29.53
##
## Marginal log-Likelihood: -42.54
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')southeast_results <- bind_rows(southeast_all %>% mutate(Region = "Southeast",Outcome = "All Injury Claims"),
southeast_cyclist %>% mutate(Region = "Southeast",Outcome = "Cyclist Injury Claims"),
southeast_pedestrian %>% mutate(Region = "Southeast",Outcome = "Pedestrian Injury Claims")
) %>%
filter(variable == "vandix_z") %>%
select(Region,Outcome,everything())Northwest
Descriptive Statistics
region | variable | n | sum | mean | sd | min | max | missing | n_zero | p_zero |
|---|---|---|---|---|---|---|---|---|---|---|
Northwest | n_claims | 78 | 4,596.0 | 58.9 | 105.9 | 0.0 | 746.0 | 0 | 5 | 6.4 |
n_casualty_claims | 498.0 | 6.4 | 11.8 | 0.0 | 67.0 | 0 | 17 | 21.8 | ||
n_cyclist_claims | 18.0 | 0.2 | 0.6 | 0.0 | 3.0 | 0 | 67 | 85.9 | ||
n_cyclist_casualty_claims | 9.0 | 0.1 | 0.4 | 0.0 | 2.0 | 0 | 70 | 89.7 | ||
n_pedestrian_claims | 51.0 | 0.7 | 1.2 | 0.0 | 5.0 | 0 | 54 | 69.2 | ||
n_pedestrian_casualty_claims | 41.0 | 0.5 | 1.1 | 0.0 | 5.0 | 0 | 57 | 73.1 | ||
population | 33,048.0 | 434.8 | 139.6 | 0.0 | 747.0 | 2 | 1 | 1.3 | ||
total_roads_km | 1,124.4 | 14.4 | 29.7 | 0.0 | 151.7 | 0 | 3 | 3.8 | ||
roads_prop_highway_arterial | 11.7 | 0.2 | 0.2 | 0.0 | 1.0 | 3 | 36 | 46.2 | ||
no_highschool_prevalance | 18.6 | 0.2 | 0.1 | 0.1 | 0.6 | 3 | 0 | 0.0 | ||
unemployment_rate | 934.5 | 12.5 | 9.3 | 0.0 | 50.0 | 3 | 4 | 5.1 | ||
hh_avg_income | 5,356,607.0 | 74,397.3 | 16,951.4 | 39,451.0 | 115,632.0 | 6 | 0 | 0.0 | ||
participation_rate | 5,081.6 | 67.8 | 7.7 | 36.4 | 81.4 | 3 | 0 | 0.0 | ||
university_degree_prevalance | 34.3 | 0.5 | 0.1 | 0.2 | 0.6 | 3 | 0 | 0.0 | ||
lone_parent_fam_prevalence | 14.7 | 0.2 | 0.1 | 0.0 | 0.4 | 3 | 0 | 0.0 | ||
home_owner_prevalence | 52.5 | 0.7 | 0.2 | 0.2 | 1.0 | 3 | 0 | 0.0 | ||
vandix | 48.8 | 0.7 | 0.9 | -0.7 | 3.5 | 6 | 0 | 0.0 | ||
ale_index | -108.7 | -1.5 | 0.5 | -2.1 | -0.5 | 5 | 0 | 0.0 | ||
canbics_index | 0.0 | 0.0 | 0.0 | 0.0 | 0.0 | 4 | 74 | 94.9 |
Maps
claims <- ggplot() +
geom_sf(data = northwest_sf, aes(fill = n_casualty_claims,colour=n_casualty_claims)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Insurance Claims",
type = "aggregation", palette = "Earth", direction = -1) +
theme_void() +
ggtitle(
"Vancouver"
)
vandix <- ggplot() +
geom_sf(data = northwest_sf, aes(fill = vandix_z_c,colour=vandix_z_c)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
scale_colour_carto_d(name = "VanDIX Score ",
type = "diverging", palette = "Earth", direction = -1) +
theme_void()
total_roads <- ggplot() +
geom_sf(data = northwest_sf, aes(fill = total_roads_km,colour=total_roads_km)) +
coord_sf(crs = "+proj=utm +zone=10 +datum=NAD83 +units=m +no_defs") +
scale_fill_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
scale_colour_carto_c(name = "Kilometres of Road",
type = "quantitative", palette = "Earth", direction = -1) +
theme_void()
cowplot::plot_grid(claims,vandix,total_roads,ncol=1)
Setting up Spatial Weights Matrix
#### Define Spatial Neighrbourhoods
northwest_sp <- as(northwest_sf, "Spatial")
northwest_sp$ui <- 1:nrow(northwest_sp@data)
coords <- coordinates(northwest_sp)
northwest_nb <- poly2nb(northwest_sp, queen = TRUE)
northwest_nb## Neighbour list object:
## Number of regions: 78
## Number of nonzero links: 406
## Percentage nonzero weights: 6.673241
## Average number of links: 5.205128#assign nearest neighbour for no links
northwest_nb <- assign_nearest_neighbors(northwest_nb,northwest_sp)
plot(northwest_sp, border = grey(0.5))
plot(northwest_nb,
coords = coords,
add = TRUE, pch = 16, lwd = 2)
Spatial Dependency in each crash type
listw <- nb2listw(northwest_nb,zero.policy = TRUE)
all_cc_mi <- moran.test(northwest_sf$n_casualty_claims, listw,zero.policy = TRUE)
all_cc_mi##
## Moran I test under randomisation
##
## data: northwest_sf$n_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 4.5914, p-value = 2.201e-06
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.290791318 -0.012987013 0.004377475cyc_cc_mi <- moran.test(northwest_sf$n_cyclist_casualty_claims, listw,zero.policy = TRUE)
cyc_cc_mi##
## Moran I test under randomisation
##
## data: northwest_sf$n_cyclist_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 1.1465, p-value = 0.1258
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.064635043 -0.012987013 0.004583381pd_cc_mi <- moran.test(northwest_sf$n_pedestrian_casualty_claims, listw,zero.policy = TRUE)
pd_cc_mi##
## Moran I test under randomisation
##
## data: northwest_sf$n_pedestrian_casualty_claims
## weights: listw
##
## Moran I statistic standard deviate = 4.4101, p-value = 5.165e-06
## alternative hypothesis: greater
## sample estimates:
## Moran I statistic Expectation Variance
## 0.298300116 -0.012987013 0.004982139BYM2 Models
nb2INLA("northwest.adj", northwest_nb)
g <- inla.read.graph(filename = "northwest.adj")
# Model Set 1: Total Casualty Claim Crashes
northwest_models_1 <- cma_models(northwest_sp@data, "Northwest", "n_casualty_claims", "vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"northwest_all <- clean_fixed_effects(northwest_models_1$fixed_effects)
map(northwest_models_1$models, ~ summary(.x))## $Unadjusted
## Time used:
## Pre = 14.1, Running = 0.384, Post = 0.152, Total = 14.6
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 0.682 0.196 0.281 0.687 1.053 0.687 0
## vandix_z 0.281 0.122 0.042 0.280 0.521 0.280 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.632 0.158 0.381 0.612 0.999 0.574
## Phi for ui 0.367 0.202 0.054 0.345 0.785 0.226
##
## Deviance Information Criterion (DIC) ...............: 354.83
## Deviance Information Criterion (DIC, saturated) ....: 151.31
## Effective number of parameters .....................: 61.55
##
## Watanabe-Akaike information criterion (WAIC) ...: 364.65
## Effective number of parameters .................: 48.83
##
## Marginal log-Likelihood: -189.93
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14, Running = 0.388, Post = 0.138, Total = 14.5
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) 0.162 0.216 -0.275 0.166 0.573 0.166 0
## vandix_z 0.362 0.115 0.138 0.361 0.591 0.361 0
## ln_roads_km_c 0.780 0.160 0.469 0.780 1.096 0.780 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 0.559 0.144 0.323 0.543 0.885 0.515
## Phi for ui 0.810 0.148 0.430 0.851 0.983 0.950
##
## Deviance Information Criterion (DIC) ...............: 338.75
## Deviance Information Criterion (DIC, saturated) ....: 135.24
## Effective number of parameters .....................: 54.97
##
## Watanabe-Akaike information criterion (WAIC) ...: 340.10
## Effective number of parameters .................: 39.26
##
## Marginal log-Likelihood: -184.44
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14, Running = 0.418, Post = 0.123, Total = 14.5
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -5.270 8.504 -21.945 -5.270 11.410
## vandix_z 0.222 0.107 0.013 0.222 0.433
## ln_roads_km_c 0.423 0.284 -0.130 0.422 0.986
## roads_prop_highway_arterial_z 0.228 0.221 -0.209 0.228 0.663
## ale_index_z 3.505 0.914 1.692 3.510 5.292
## canbics_index_z -10.849 9.886 -30.222 -10.854 8.552
## population_100_c 0.372 0.123 0.130 0.373 0.614
## mode kld
## (Intercept) -5.270 0
## vandix_z 0.222 0
## ln_roads_km_c 0.422 0
## roads_prop_highway_arterial_z 0.228 0
## ale_index_z 3.510 0
## canbics_index_z -10.854 0
## population_100_c 0.373 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1.159 0.296 0.688 1.121 1.846 1.049
## Phi for ui 0.192 0.214 0.002 0.104 0.776 0.002
##
## Deviance Information Criterion (DIC) ...............: 353.41
## Deviance Information Criterion (DIC, saturated) ....: 149.90
## Effective number of parameters .....................: 61.24
##
## Watanabe-Akaike information criterion (WAIC) ...: 366.70
## Effective number of parameters .................: 50.50
##
## Marginal log-Likelihood: -186.55
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 2: Total Casualty Cyclist Claim Crashes
northwest_models_2 <- cma_models(northwest_sp@data,"Northwest","n_cyclist_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"northwest_cyclist <- clean_fixed_effects(northwest_models_2$fixed_effects)
map(northwest_models_2$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 13.9, Running = 0.419, Post = 0.133, Total = 14.5
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -2.583 0.466 -3.500 -2.582 -1.670 -2.582 0
## vandix_z 0.237 0.210 -0.173 0.237 0.649 0.237 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1271.475 14399.03 2.365 86.105 7984.694 3.618
## Phi for ui 0.343 0.28 0.009 0.269 0.927 0.012
##
## Deviance Information Criterion (DIC) ...............: 60.72
## Deviance Information Criterion (DIC, saturated) ....: 44.01
## Effective number of parameters .....................: 2.28
##
## Watanabe-Akaike information criterion (WAIC) ...: 60.86
## Effective number of parameters .................: 2.37
##
## Marginal log-Likelihood: 6.05
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 14, Running = 0.368, Post = 0.142, Total = 14.5
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -2.630 0.519 -3.651 -2.630 -1.615 -2.630 0
## vandix_z 0.246 0.217 -0.180 0.246 0.672 0.246 0
## ln_roads_km_c -0.063 0.358 -0.765 -0.063 0.638 -0.063 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1284.692 14712.83 2.325 85.414 8040.646 3.540
## Phi for ui 0.343 0.28 0.009 0.269 0.927 0.012
##
## Deviance Information Criterion (DIC) ...............: 62.65
## Deviance Information Criterion (DIC, saturated) ....: 45.93
## Effective number of parameters .....................: 3.17
##
## Watanabe-Akaike information criterion (WAIC) ...: 62.90
## Effective number of parameters .................: 3.26
##
## Marginal log-Likelihood: 1.57
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 14, Running = 0.403, Post = 0.124, Total = 14.5
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -56.307 30.490 -97.521 -65.611 6.254
## vandix_z 0.069 0.370 -0.498 0.030 1.253
## ln_roads_km_c -0.993 1.583 -5.483 -0.811 2.554
## roads_prop_highway_arterial_z 0.654 0.868 -1.945 0.723 2.010
## ale_index_z 2.985 3.395 -7.643 3.333 8.266
## canbics_index_z -63.899 35.644 -113.806 -73.758 6.584
## population_100_c 0.758 0.628 -0.010 0.606 2.551
## mode kld
## (Intercept) -71.155 0.000
## vandix_z -0.003 0.005
## ln_roads_km_c -0.881 0.024
## roads_prop_highway_arterial_z 0.816 0.001
## ale_index_z 3.650 0.003
## canbics_index_z -88.299 0.000
## population_100_c 0.579 0.004
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1381.138 16139.07 2.442 88.177 8570.257 3.794
## Phi for ui 0.344 0.28 0.009 0.269 0.927 0.012
##
## Deviance Information Criterion (DIC) ...............: -1.29e+13
## Deviance Information Criterion (DIC, saturated) ....: -1.29e+13
## Effective number of parameters .....................: -1.29e+13
##
## Watanabe-Akaike information criterion (WAIC) ...: 1835.70
## Effective number of parameters .................: 889.55
##
## Marginal log-Likelihood: -7.91
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')# Model Set 3: Total Casualty Cyclist Claim Crashes
northwest_models_3 <- cma_models(northwest_sp@data,"Northwest","n_pedestrian_casualty_claims","vandix_z") ## [1] "Fit unadjusted"
## [1] "Fit adjusted by road length"
## [1] "Fit adjusted by road length and covariates"northwest_pedestrian <- clean_fixed_effects(northwest_models_3$fixed_effects)
map(northwest_models_3$models,~summary(.x))## $Unadjusted
## Time used:
## Pre = 14.3, Running = 0.387, Post = 0.118, Total = 14.8
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.469 0.344 -2.206 -1.450 -0.848 -1.400 0
## vandix_z 0.275 0.144 -0.006 0.274 0.563 0.274 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1.544 0.898 0.516 1.324 3.903 0.993
## Phi for ui 0.341 0.234 0.027 0.295 0.849 0.080
##
## Deviance Information Criterion (DIC) ...............: 142.59
## Deviance Information Criterion (DIC, saturated) ....: 90.73
## Effective number of parameters .....................: 19.88
##
## Watanabe-Akaike information criterion (WAIC) ...: 147.89
## Effective number of parameters .................: 19.43
##
## Marginal log-Likelihood: -45.03
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted1
## Time used:
## Pre = 13.9, Running = 0.382, Post = 0.137, Total = 14.4
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant mode kld
## (Intercept) -1.681 0.402 -2.536 -1.661 -0.951 -1.602 0
## vandix_z 0.303 0.151 0.010 0.302 0.606 0.302 0
## ln_roads_km_c 0.252 0.242 -0.212 0.247 0.740 0.247 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 1.37 0.782 0.461 1.179 3.417 0.892
## Phi for ui 0.45 0.253 0.050 0.433 0.917 0.205
##
## Deviance Information Criterion (DIC) ...............: 141.37
## Deviance Information Criterion (DIC, saturated) ....: 89.50
## Effective number of parameters .....................: 21.00
##
## Watanabe-Akaike information criterion (WAIC) ...: 145.97
## Effective number of parameters .................: 19.60
##
## Marginal log-Likelihood: -49.40
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
##
##
## $Adjusted2
## Time used:
## Pre = 13.9, Running = 0.405, Post = 0.128, Total = 14.4
## Fixed effects:
## mean sd 0.025quant 0.5quant 0.975quant
## (Intercept) -5.513 9.156 -23.468 -5.513 12.442
## vandix_z 0.175 0.153 -0.126 0.174 0.478
## ln_roads_km_c 0.143 0.441 -0.699 0.134 1.037
## roads_prop_highway_arterial_z 0.258 0.354 -0.452 0.263 0.941
## ale_index_z 4.087 1.278 1.588 4.078 6.636
## canbics_index_z -9.180 10.658 -30.075 -9.181 11.724
## population_100_c 0.349 0.162 0.032 0.347 0.674
## mode kld
## (Intercept) -5.513 0
## vandix_z 0.174 0
## ln_roads_km_c 0.134 0
## roads_prop_highway_arterial_z 0.262 0
## ale_index_z 4.078 0
## canbics_index_z -9.181 0
## population_100_c 0.347 0
##
## Random effects:
## Name Model
## ui BYM2 model
##
## Model hyperparameters:
## mean sd 0.025quant 0.5quant 0.975quant mode
## Precision for ui 2.332 1.734 0.617 1.858 6.956 1.247
## Phi for ui 0.305 0.218 0.024 0.255 0.804 0.071
##
## Deviance Information Criterion (DIC) ...............: 159.78
## Deviance Information Criterion (DIC, saturated) ....: 107.92
## Effective number of parameters .....................: 28.84
##
## Watanabe-Akaike information criterion (WAIC) ...: 189.93
## Effective number of parameters .................: 39.10
##
## Marginal log-Likelihood: -55.76
## CPO, PIT is computed
## Posterior summaries for the linear predictor and the fitted values are computed
## (Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')northwest_results <- bind_rows(northwest_all %>% mutate(Region = "Northwest",Outcome = "All Injury Claims"),
northwest_cyclist %>% mutate(Region = "Northwest",Outcome = "Cyclist Injury Claims"),
northwest_pedestrian %>% mutate(Region = "Northwest",Outcome = "Pedestrian Injury Claims")
) %>%
filter(variable == "vandix_z") %>%
select(Region,Outcome,everything())
beepr::beep(8)Results
Study Area
## Downloading: 12 kB Downloading: 12 kB Downloading: 12 kB Downloading: 12 kB Downloading: 16 kB Downloading: 16 kB Downloading: 16 kB Downloading: 16 kB Downloading: 24 kB Downloading: 24 kB Downloading: 24 kB Downloading: 24 kB
Number of DAs | Total Injury Claims | Injury Claims Mean (SD) | Total Cycling Injury Claims | Cycling Injury Claims Mean (SD) | Total Pedestrian Injury Claims | Pedestrian Injury Claims Mean (SD) | Population | VANDIX Mean (SD) | Total Roads (km) | Highway/Arterial Roads Proportion Mean (SD) | No Highschool Prevalence Mean (SD) | Unemployment Rate Mean (SD) | Average Household Income Mean (SD) | Participation Rate Mean (SD) | University Degree Prevalence Mean (SD) | Lone Parent Family Prevalence Mean (SD) | Home Ownership Prevalence Mean (SD) | ALE Index Mean (SD) | CANBICS Index Mean (SD) |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
6,503 | 211,852 | 32.6 (70.8) | 7,287 | 1.12 (2.77) | 10,531 | 1.62 (3.46) | 4,477,748 | -0.016 (0.64) | 47,865.21 | 0.146 (0.174) | 0.158 (0.0805) | 6.74 (4.62) | 80100 (30400) | 63.8 (10.9) | 0.546 (0.118) | 0.16 (0.0787) | 0.691 (0.222) | 0.692 (3.07) | 3.15 (3.65) |
The analysis included 6503 dissemination areas across 20 different metropolitan regions with a total of 211852 motor-vehicle crashes resulting in at least one injury, and 7287 such crashes involved at least one cyclist, while 10531 involved at least one pedestrian.
region_name | Number of DAs | Total Injury Claims | Injury Claims Mean (SD) | Total Cycling Injury Claims | Cycling Injury Claims Mean (SD) | Total Pedestrian Injury Claims | Pedestrian Injury Claims Mean (SD) | Population | VANDIX Mean (SD) | Total Roads (km) | Highway/Arterial Roads Proportion Mean (SD) | No Highschool Prevalence Mean (SD) | Unemployment Rate Mean (SD) | Average Household Income Mean (SD) | Participation Rate Mean (SD) | University Degree Prevalence Mean (SD) | Lone Parent Family Prevalence Mean (SD) | Home Ownership Prevalence Mean (SD) | ALE Index Mean (SD) | CANBICS Index Mean (SD) |
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
Vancouver-Squamish | 3,622 | 149,692 | 41.3 (84.8) | 4,968 | 1.37 (3.24) | 7,502 | 2.07 (4.02) | 2,667,057 | -0.14 (0.58) | 15,056.825 | 0.159 (0.184) | 0.143 (0.0767) | 5.88 (3.58) | 85400 (34200) | 64.9 (10.1) | 0.569 (0.121) | 0.155 (0.0659) | 0.671 (0.23) | 1.63 (3.58) | 4.57 (4.11) |
Victoria | 581 | 12,338 | 21.2 (34.2) | 1,027 | 1.77 (3.17) | 715 | 1.23 (2.54) | 397,237 | -0.24 (0.52) | 3,261.164 | 0.116 (0.14) | 0.121 (0.0581) | 5.63 (3.87) | 77300 (26900) | 63.6 (11.1) | 0.597 (0.0978) | 0.153 (0.08) | 0.654 (0.231) | 0.607 (1.78) | 2.12 (1.71) |
Central Island-Powell River | 611 | 10,762 | 17.6 (28.3) | 271 | 0.444 (0.938) | 547 | 0.895 (1.96) | 357,193 | 0.2 (0.65) | 5,652.122 | 0.13 (0.161) | 0.177 (0.0756) | 8.35 (5.15) | 65700 (17200) | 56.5 (11) | 0.524 (0.0938) | 0.167 (0.0948) | 0.744 (0.188) | -1.07 (0.719) | 0.588 (0.914) |
Okanagan | 466 | 13,236 | 28.4 (66.6) | 485 | 1.04 (2.34) | 564 | 1.21 (3.23) | 336,628 | 0.12 (0.6) | 4,759.053 | 0.13 (0.162) | 0.164 (0.0651) | 7.76 (4.73) | 73000 (26800) | 61.3 (12.9) | 0.527 (0.0836) | 0.158 (0.0831) | 0.728 (0.208) | -0.843 (0.959) | 1.9 (1.81) |
Fraser Valley | 472 | 15,730 | 33.3 (65.2) | 360 | 0.763 (1.58) | 725 | 1.54 (2.96) | 309,493 | 0.31 (0.6) | 3,678.908 | 0.0858 (0.138) | 0.216 (0.0773) | 7.05 (5.06) | 75400 (21100) | 63.6 (11.6) | 0.444 (0.088) | 0.166 (0.101) | 0.741 (0.172) | -0.24 (1.73) | 1.39 (1.3) |
North Central | 359 | 4,490 | 12.5 (25.2) | 54 | 0.15 (0.544) | 203 | 0.565 (1.62) | 182,818 | 0.5 (0.71) | 9,989.620 | 0.13 (0.145) | 0.227 (0.0863) | 10.5 (6.84) | 77700 (21900) | 67.9 (9.34) | 0.453 (0.0895) | 0.182 (0.1) | 0.711 (0.236) | -1.32 (0.637) | 0.582 (0.981) |
Kamloops-Salmon Arm | 204 | 4,069 | 19.9 (39.5) | 89 | 0.436 (1.06) | 183 | 0.897 (2.28) | 133,847 | 0.18 (0.67) | 2,744.521 | 0.245 (0.19) | 0.17 (0.0717) | 8.32 (5.27) | 71500 (20400) | 62.3 (9.95) | 0.511 (0.0944) | 0.178 (0.101) | 0.74 (0.201) | -0.714 (1.22) | 1.4 (1.47) |
Southeast | 110 | 1,037 | 9.43 (17.8) | 24 | 0.218 (0.612) | 51 | 0.464 (1.23) | 60,427 | 0.079 (0.58) | 1,598.577 | 0.128 (0.156) | 0.156 (0.0625) | 8.06 (4.35) | 68700 (16800) | 62.2 (8.54) | 0.549 (0.0921) | 0.165 (0.0835) | 0.743 (0.203) | -1.24 (0.743) | 0.321 (0.535) |
Northwest | 78 | 498 | 6.38 (11.8) | 9 | 0.115 (0.36) | 41 | 0.526 (1.05) | 33,048 | 0.68 (0.91) | 1,124.422 | 0.156 (0.223) | 0.247 (0.098) | 12.5 (9.31) | 74400 (17000) | 67.8 (7.67) | 0.457 (0.102) | 0.196 (0.0949) | 0.7 (0.189) | -1.49 (0.483) | 0 (0) |
Association between Deprivation Index and Traffic Injury
Crude Rates
All Motor Vehicle Traffic Injury Crashes
Bicycle-Motor Vehicle Traffic Injury Crashes
Pedestrian-Motor Vehicle Traffic Injury Crashes
Spatial Models
All Motor Vehicle Traffic Injury Crashes
Variable | Region | Unadjusted IRR (95% CI) | Fully Adjusted IRR (95% CI)a |
|---|---|---|---|
VanDIX (SD) | Vancouver-Squamish | 1.24 (1.16, 1.32) | 1.25 (1.2, 1.3) |
Victoria | 1.22 (1.07, 1.38) | 1.22 (1.11, 1.35) | |
Central Island-Powell River | 1.34 (1.19, 1.52) | 1.37 (1.26, 1.5) | |
Okanagan | 1.41 (1.18, 1.69) | 1.43 (1.26, 1.62) | |
Fraser Valley | 1.26 (1.07, 1.49) | 1.15 (1.02, 1.29) | |
North Central | 1.21 (1.04, 1.41) | 1.22 (1.07, 1.39) | |
Kamloops-Salmon Arm | 1.28 (1.02, 1.61) | 1.23 (1.04, 1.46) | |
Southeast | 1.42 (1.06, 1.91) | 1.42 (1.1, 1.84) | |
Northwest | 1.32 (1.04, 1.68) | 1.25 (1.01, 1.54) | |
aAdjusted for (i) total length of roads; (ii) proportion of roads classified as highway or arterial; (iii) the Canadian Active Living Environment Index; (iv) the Canadian Bikeway Comfort and Safety Index; (v) Total population | |||
Bicycle-Motor Vehicle Traffic Injury Crashes
Variable | Region | Unadjusted IRR (95% CI) | Fully Adjusted IRR (95% CI)a |
|---|---|---|---|
VanDIX (SD) | Vancouver-Squamish | 1.16 (1.08, 1.25) | 1.2 (1.13, 1.28) |
Victoria | 1.16 (0.983, 1.37) | 1.13 (0.962, 1.32) | |
Central Island-Powell River | 1.4 (1.21, 1.62) | 1.42 (1.22, 1.66) | |
Okanagan | 1.9 (1.54, 2.36) | 1.76 (1.45, 2.15) | |
Fraser Valley | 1.4 (1.12, 1.75) | 1.08 (0.87, 1.35) | |
North Central | 1.36 (1.04, 1.79) | 1.43 (1.07, 1.92) | |
Kamloops-Salmon Arm | 1.45 (1.07, 1.95) | 1.48 (1.03, 2.07) | |
Southeast | 1.82 (1.17, 2.89) | 2.29 (1.4, 3.76) | |
Northwest | 1.27 (0.841, 1.91) | 1.07 (0.607, 3.5) | |
aAdjusted for (i) total length of roads; (ii) proportion of roads classified as highway or arterial; (iii) the Canadian Active Living Environment Index; (iv) the Canadian Bikeway Comfort and Safety Index; (v) Total population | |||
Pedestrian-Motor Vehicle Traffic Injury Crashes
Variable | Region | Unadjusted IRR (95% CI) | Fully Adjusted IRR (95% CI)a |
|---|---|---|---|
VanDIX (SD) | Vancouver-Squamish | 1.3 (1.21, 1.4) | 1.31 (1.23, 1.39) |
Victoria | 1.42 (1.2, 1.67) | 1.37 (1.17, 1.61) | |
Central Island-Powell River | 1.74 (1.5, 2.03) | 1.68 (1.45, 1.95) | |
Okanagan | 1.95 (1.57, 2.42) | 1.75 (1.4, 2.18) | |
Fraser Valley | 1.65 (1.37, 1.99) | 1.32 (1.11, 1.57) | |
North Central | 1.47 (1.19, 1.8) | 1.4 (1.13, 1.73) | |
Kamloops-Salmon Arm | 1.92 (1.41, 2.64) | 1.94 (1.44, 2.63) | |
Southeast | 1.72 (1.08, 2.78) | 1.83 (1.12, 3.02) | |
Northwest | 1.32 (0.994, 1.76) | 1.19 (0.882, 1.61) | |
aAdjusted for (i) total length of roads; (ii) proportion of roads classified as highway or arterial; (iii) the Canadian Active Living Environment Index; (iv) the Canadian Bikeway Comfort and Safety Index; (v) Total population | |||
Supplementary Material
Built Environment by VanDIX
- Of the 9 regions, we observed 5 had on average a positive relationship with deprivation index and the proportion of roads in the neighborhood that were arterial or highway
- Of the 9 regions, we observed 8 had on average a positive relationship with deprivation index and thactive living environment index score
- Of the 9 regions, we observed 8 had on average a positive relationship with deprivation index and increased access to better bicycle facilities.
Model Details
Vancouver Model Details
All
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | 2.800 (2.770, 2.830) | 3.250 (3.220, 3.280) | 2.990 (2.960, 3.020) |
vandix_z | 0.211 (0.147, 0.276) | 0.327 (0.276, 0.378) | 0.222 (0.180, 0.265) |
ln_roads_km_c | 1.290 (1.240, 1.350) | 0.990 (0.937, 1.040) | |
roads_prop_highway_arterial_z |
| 0.593 (0.562, 0.624) | |
ale_index_z |
| 0.109 (0.036, 0.182) | |
canbics_index_z |
| 0.209 (0.137, 0.282) | |
population_100_c |
| 0.017 (0.012, 0.022) | |
Hyper Parameters | |||
Precision for ui | 0.311 (0.284, 0.339) | 0.508 (0.470, 0.548) | 0.838 (0.765, 0.915) |
Phi for ui | 0.794 (0.751, 0.833) | 0.811 (0.778, 0.841) | 0.777 (0.733, 0.818) |
Cyclists
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -0.523 (-0.588, -0.461) | -0.134 (-0.186, -0.083) | -0.291 (-0.352, -0.231) |
vandix_z | 0.148 (0.075, 0.221) | 0.217 (0.154, 0.280) | 0.185 (0.123, 0.248) |
ln_roads_km_c | 1.070 (1.000, 1.140) | 0.844 (0.767, 0.921) | |
roads_prop_highway_arterial_z |
| 0.301 (0.254, 0.348) | |
ale_index_z |
| 0.041 (-0.048, 0.130) | |
canbics_index_z |
| 0.092 (-0.001, 0.186) | |
population_100_c |
| 0.022 (0.016, 0.028) | |
Hyper Parameters | |||
Precision for ui | 0.567 (0.489, 0.653) | 0.622 (0.532, 0.718) | 0.822 (0.692, 0.964) |
Phi for ui | 0.633 (0.534, 0.724) | 0.934 (0.874, 0.974) | 0.871 (0.790, 0.931) |
Pedestrians
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -0.035 (-0.089, 0.018) | 0.319 (0.275, 0.363) | 0.115 (0.060, 0.168) |
vandix_z | 0.262 (0.191, 0.333) | 0.304 (0.242, 0.367) | 0.268 (0.208, 0.328) |
ln_roads_km_c | 0.997 (0.930, 1.060) | 0.717 (0.641, 0.793) | |
roads_prop_highway_arterial_z |
| 0.394 (0.351, 0.438) | |
ale_index_z |
| 0.232 (0.141, 0.322) | |
canbics_index_z |
| 0.042 (-0.051, 0.135) | |
population_100_c |
| 0.027 (0.021, 0.034) | |
Hyper Parameters | |||
Precision for ui | 0.496 (0.430, 0.567) | 0.450 (0.396, 0.506) | 0.713 (0.602, 0.835) |
Phi for ui | 0.663 (0.573, 0.746) | 0.950 (0.907, 0.979) | 0.850 (0.763, 0.917) |
Victoria Model Details
All
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | 2.340 (2.270, 2.410) | 2.360 (2.310, 2.410) | 2.510 (2.420, 2.600) |
vandix_z | 0.196 (0.068, 0.325) | 0.326 (0.222, 0.430) | 0.202 (0.104, 0.301) |
ln_roads_km_c | 1.240 (1.100, 1.380) | 1.010 (0.850, 1.160) | |
roads_prop_highway_arterial_z |
| 0.462 (0.367, 0.558) | |
ale_index_z |
| 0.409 (0.081, 0.739) | |
canbics_index_z |
| 0.340 (0.032, 0.646) | |
population_100_c |
| 0.015 (-0.002, 0.031) | |
Hyper Parameters | |||
Precision for ui | 0.490 (0.390, 0.599) | 0.642 (0.547, 0.744) | 0.845 (0.679, 1.020) |
Phi for ui | 0.857 (0.721, 0.947) | 0.981 (0.930, 0.999) | 0.948 (0.838, 0.994) |
Cyclists
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -0.141 (-0.282, -0.006) | -0.144 (-0.273, -0.018) | 0.039 (-0.120, 0.193) |
vandix_z | 0.150 (-0.017, 0.318) | 0.250 (0.094, 0.406) | 0.118 (-0.038, 0.275) |
ln_roads_km_c | 1.100 (0.871, 1.320) | 0.887 (0.635, 1.140) | |
roads_prop_highway_arterial_z |
| 0.308 (0.165, 0.451) | |
ale_index_z |
| 0.401 (-0.033, 0.838) | |
canbics_index_z |
| 0.587 (0.169, 1.010) | |
population_100_c |
| 0.030 (0.004, 0.058) | |
Hyper Parameters | |||
Precision for ui | 0.679 (0.511, 0.882) | 0.663 (0.489, 0.868) | 0.867 (0.622, 1.160) |
Phi for ui | 0.623 (0.414, 0.806) | 0.902 (0.739, 0.984) | 0.828 (0.616, 0.957) |
Pedestrians
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -0.439 (-0.589, -0.296) | -0.430 (-0.573, -0.292) | -0.313 (-0.493, -0.141) |
vandix_z | 0.350 (0.184, 0.516) | 0.452 (0.292, 0.612) | 0.317 (0.156, 0.478) |
ln_roads_km_c | 0.854 (0.625, 1.090) | 0.578 (0.313, 0.845) | |
roads_prop_highway_arterial_z |
| 0.384 (0.228, 0.539) | |
ale_index_z |
| 0.540 (0.080, 1.010) | |
canbics_index_z |
| 0.326 (-0.136, 0.786) | |
population_100_c |
| 0.037 (0.010, 0.064) | |
Hyper Parameters | |||
Precision for ui | 0.633 (0.458, 0.847) | 0.619 (0.458, 0.816) | 0.758 (0.543, 1.020) |
Phi for ui | 0.791 (0.569, 0.937) | 0.950 (0.822, 0.996) | 0.924 (0.745, 0.993) |
Central Island - Powell River Model Details
All
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | 1.930 (1.830, 2.020) | 1.360 (1.240, 1.470) | 1.960 (1.730, 2.190) |
vandix_z | 0.295 (0.175, 0.416) | 0.435 (0.330, 0.539) | 0.317 (0.227, 0.408) |
ln_roads_km_c | 1.020 (0.888, 1.160) | 0.650 (0.506, 0.794) | |
roads_prop_highway_arterial_z |
| 0.668 (0.575, 0.761) | |
ale_index_z |
| 1.420 (0.937, 1.900) | |
canbics_index_z |
| -0.780 (-1.230, -0.328) | |
population_100_c |
| 0.063 (0.037, 0.090) | |
Hyper Parameters | |||
Precision for ui | 0.410 (0.309, 0.537) | 0.442 (0.345, 0.558) | 0.670 (0.517, 0.852) |
Phi for ui | 0.571 (0.383, 0.741) | 0.736 (0.613, 0.836) | 0.726 (0.598, 0.832) |
Cyclists
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -1.430 (-1.650, -1.220) | -1.730 (-2.010, -1.460) | -1.150 (-1.570, -0.753) |
vandix_z | 0.335 (0.188, 0.483) | 0.399 (0.250, 0.548) | 0.354 (0.203, 0.504) |
ln_roads_km_c | 0.506 (0.280, 0.737) | 0.292 (0.018, 0.568) | |
roads_prop_highway_arterial_z |
| 0.431 (0.272, 0.592) | |
ale_index_z |
| 1.160 (0.350, 1.970) | |
canbics_index_z |
| -0.321 (-1.040, 0.394) | |
population_100_c |
| 0.060 (0.014, 0.107) | |
Hyper Parameters | |||
Precision for ui | 1.060 (0.680, 1.600) | 0.927 (0.561, 1.460) | 1.180 (0.691, 1.910) |
Phi for ui | 0.320 (0.089, 0.626) | 0.551 (0.249, 0.823) | 0.496 (0.197, 0.788) |
Pedestrians
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -1.250 (-1.490, -1.040) | -1.600 (-1.880, -1.330) | -0.755 (-1.170, -0.361) |
vandix_z | 0.556 (0.405, 0.710) | 0.611 (0.456, 0.769) | 0.518 (0.370, 0.666) |
ln_roads_km_c | 0.653 (0.401, 0.908) | 0.328 (0.052, 0.605) | |
roads_prop_highway_arterial_z |
| 0.552 (0.394, 0.713) | |
ale_index_z |
| 2.370 (1.610, 3.150) | |
canbics_index_z |
| -1.150 (-1.840, -0.468) | |
population_100_c |
| 0.115 (0.071, 0.160) | |
Hyper Parameters | |||
Precision for ui | 0.605 (0.439, 0.830) | 0.468 (0.302, 0.701) | 0.749 (0.489, 1.120) |
Phi for ui | 0.206 (0.030, 0.511) | 0.587 (0.295, 0.830) | 0.446 (0.166, 0.736) |
Okanagan Model Details
All
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | 2.160 (2.070, 2.260) | 1.580 (1.480, 1.680) | 1.960 (1.800, 2.130) |
vandix_z | 0.344 (0.163, 0.525) | 0.464 (0.318, 0.610) | 0.355 (0.228, 0.484) |
ln_roads_km_c | 1.150 (1.000, 1.290) | 0.949 (0.801, 1.100) | |
roads_prop_highway_arterial_z |
| 0.605 (0.499, 0.711) | |
ale_index_z |
| 0.510 (-0.006, 1.030) | |
canbics_index_z |
| 0.016 (-0.406, 0.438) | |
population_100_c |
| 0.031 (0.014, 0.048) | |
Hyper Parameters | |||
Precision for ui | 0.237 (0.172, 0.317) | 0.241 (0.200, 0.286) | 0.339 (0.271, 0.412) |
Phi for ui | 0.843 (0.704, 0.934) | 0.985 (0.950, 0.998) | 0.979 (0.929, 0.998) |
Cyclists
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -1.110 (-1.340, -0.903) | -1.510 (-1.790, -1.240) | -1.010 (-1.350, -0.697) |
vandix_z | 0.642 (0.429, 0.858) | 0.705 (0.500, 0.913) | 0.566 (0.369, 0.765) |
ln_roads_km_c | 0.721 (0.484, 0.961) | 0.616 (0.348, 0.887) | |
roads_prop_highway_arterial_z |
| 0.395 (0.221, 0.569) | |
ale_index_z |
| 0.703 (-0.040, 1.450) | |
canbics_index_z |
| 0.354 (-0.220, 0.925) | |
population_100_c |
| 0.036 (0.008, 0.065) | |
Hyper Parameters | |||
Precision for ui | 0.426 (0.269, 0.634) | 0.339 (0.230, 0.474) | 0.428 (0.288, 0.608) |
Phi for ui | 0.788 (0.568, 0.933) | 0.951 (0.840, 0.995) | 0.961 (0.858, 0.997) |
Pedestrians
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -1.040 (-1.280, -0.822) | -1.440 (-1.720, -1.170) | -0.861 (-1.180, -0.554) |
vandix_z | 0.668 (0.452, 0.884) | 0.731 (0.497, 0.965) | 0.558 (0.338, 0.780) |
ln_roads_km_c | 0.805 (0.541, 1.070) | 0.632 (0.337, 0.927) | |
roads_prop_highway_arterial_z |
| 0.575 (0.404, 0.748) | |
ale_index_z |
| 1.290 (0.531, 2.070) | |
canbics_index_z |
| -0.245 (-0.877, 0.384) | |
population_100_c |
| 0.049 (0.019, 0.081) | |
Hyper Parameters | |||
Precision for ui | 0.525 (0.374, 0.720) | 0.372 (0.232, 0.561) | 0.508 (0.310, 0.786) |
Phi for ui | 0.328 (0.121, 0.593) | 0.751 (0.489, 0.925) | 0.707 (0.408, 0.911) |
Fraser Valley Model Details
All
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | 2.480 (2.400, 2.560) | 2.270 (2.180, 2.350) | 2.720 (2.600, 2.840) |
vandix_z | 0.232 (0.065, 0.399) | 0.285 (0.156, 0.415) | 0.140 (0.022, 0.257) |
ln_roads_km_c | 1.100 (0.972, 1.240) | 0.862 (0.721, 1.000) | |
roads_prop_highway_arterial_z |
| 0.580 (0.464, 0.695) | |
ale_index_z |
| 0.245 (-0.004, 0.494) | |
canbics_index_z |
| 0.183 (-0.097, 0.458) | |
population_100_c |
| 0.027 (0.003, 0.051) | |
Hyper Parameters | |||
Precision for ui | 0.299 (0.258, 0.344) | 0.712 (0.585, 0.856) | 0.953 (0.760, 1.170) |
Phi for ui | 0.767 (0.001, 1.000) | 0.802 (0.680, 0.893) | 0.805 (0.658, 0.913) |
Cyclists
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -1.220 (-1.480, -0.984) | -1.360 (-1.610, -1.120) | -0.791 (-1.080, -0.513) |
vandix_z | 0.338 (0.113, 0.557) | 0.286 (0.063, 0.509) | 0.081 (-0.139, 0.299) |
ln_roads_km_c | 0.753 (0.505, 1.000) | 0.467 (0.207, 0.731) | |
roads_prop_highway_arterial_z |
| 0.601 (0.404, 0.800) | |
ale_index_z |
| 0.695 (0.276, 1.110) | |
canbics_index_z |
| 0.091 (-0.494, 0.678) | |
population_100_c |
| 0.046 (0.001, 0.091) | |
Hyper Parameters | |||
Precision for ui | 0.679 (0.485, 0.930) | 0.623 (0.427, 0.874) | 0.897 (0.588, 1.320) |
Phi for ui | 0.315 (0.098, 0.593) | 0.752 (0.479, 0.934) | 0.611 (0.278, 0.882) |
Pedestrians
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -0.585 (-0.779, -0.402) | -0.730 (-0.915, -0.552) | -0.210 (-0.437, 0.010) |
vandix_z | 0.503 (0.318, 0.689) | 0.447 (0.269, 0.627) | 0.279 (0.109, 0.450) |
ln_roads_km_c | 0.813 (0.608, 1.020) | 0.421 (0.199, 0.646) | |
roads_prop_highway_arterial_z |
| 0.619 (0.459, 0.781) | |
ale_index_z |
| 0.488 (0.140, 0.837) | |
canbics_index_z |
| 0.025 (-0.445, 0.494) | |
population_100_c |
| 0.059 (0.023, 0.097) | |
Hyper Parameters | |||
Precision for ui | 0.740 (0.558, 0.965) | 0.660 (0.476, 0.882) | 1.040 (0.732, 1.440) |
Phi for ui | 0.372 (0.166, 0.610) | 0.852 (0.638, 0.970) | 0.675 (0.392, 0.892) |
North Central Model Details
All
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | 1.320 (1.160, 1.480) | 0.531 (0.261, 0.795) | 0.969 (0.535, 1.400) |
vandix_z | 0.191 (0.042, 0.340) | 0.272 (0.132, 0.413) | 0.201 (0.071, 0.332) |
ln_roads_km_c | 0.698 (0.521, 0.876) | 0.542 (0.371, 0.713) | |
roads_prop_highway_arterial_z |
| 0.616 (0.447, 0.784) | |
ale_index_z |
| -0.263 (-0.837, 0.313) | |
canbics_index_z |
| 0.384 (-0.228, 0.994) | |
population_100_c |
| 0.034 (-0.008, 0.076) | |
Hyper Parameters | |||
Precision for ui | 0.296 (0.218, 0.391) | 0.229 (0.172, 0.294) | 0.269 (0.204, 0.341) |
Phi for ui | 0.788 (0.650, 0.890) | 0.950 (0.871, 0.990) | 0.962 (0.890, 0.994) |
Cyclists
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -3.130 (-3.910, -2.500) | -3.340 (-4.190, -2.590) | -2.520 (-3.540, -1.580) |
vandix_z | 0.311 (0.043, 0.584) | 0.350 (0.069, 0.637) | 0.356 (0.068, 0.651) |
ln_roads_km_c | 0.160 (-0.154, 0.483) | 0.211 (-0.141, 0.570) | |
roads_prop_highway_arterial_z |
| 0.504 (0.115, 0.898) | |
ale_index_z |
| -0.474 (-1.870, 0.929) | |
canbics_index_z |
| 2.030 (0.911, 3.170) | |
population_100_c |
| -0.047 (-0.197, 0.104) | |
Hyper Parameters | |||
Precision for ui | 0.649 (0.347, 1.130) | 0.637 (0.340, 1.110) | 0.833 (0.399, 1.590) |
Phi for ui | 0.158 (0.024, 0.438) | 0.199 (0.032, 0.522) | 0.306 (0.063, 0.682) |
Pedestrians
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -1.930 (-2.340, -1.570) | -2.060 (-2.570, -1.600) | -2.110 (-3.010, -1.270) |
vandix_z | 0.382 (0.177, 0.589) | 0.392 (0.182, 0.604) | 0.335 (0.126, 0.548) |
ln_roads_km_c | 0.120 (-0.147, 0.409) | -0.063 (-0.372, 0.265) | |
roads_prop_highway_arterial_z |
| 0.650 (0.364, 0.939) | |
ale_index_z |
| -1.240 (-2.500, -0.021) | |
canbics_index_z |
| 0.642 (-0.340, 1.640) | |
population_100_c |
| 0.022 (-0.064, 0.109) | |
Hyper Parameters | |||
Precision for ui | 0.548 (0.346, 0.833) | 0.512 (0.302, 0.818) | 0.522 (0.291, 0.864) |
Phi for ui | 0.359 (0.106, 0.681) | 0.451 (0.131, 0.801) | 0.560 (0.204, 0.875) |
Kamloops-Salmon Arm
All
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | 1.950 (1.790, 2.110) | 1.430 (1.230, 1.620) | 1.130 (0.898, 1.370) |
vandix_z | 0.250 (0.024, 0.476) | 0.258 (0.055, 0.461) | 0.208 (0.042, 0.375) |
ln_roads_km_c | 0.850 (0.644, 1.060) | 0.491 (0.302, 0.681) | |
roads_prop_highway_arterial_z |
| 0.534 (0.424, 0.645) | |
ale_index_z |
| 0.091 (-0.413, 0.594) | |
canbics_index_z |
| -0.712 (-1.200, -0.224) | |
population_100_c |
| 0.099 (0.062, 0.136) | |
Hyper Parameters | |||
Precision for ui | 0.647 (0.498, 0.825) | 0.776 (0.590, 0.999) | 1.180 (0.871, 1.540) |
Phi for ui | 0.421 (0.223, 0.639) | 0.625 (0.425, 0.802) | 0.819 (0.620, 0.948) |
Cyclists
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -1.620 (-2.050, -1.230) | -2.040 (-2.620, -1.520) | -1.920 (-2.590, -1.320) |
vandix_z | 0.372 (0.064, 0.668) | 0.394 (0.069, 0.709) | 0.391 (0.034, 0.726) |
ln_roads_km_c | 0.586 (0.200, 0.990) | 0.195 (-0.224, 0.624) | |
roads_prop_highway_arterial_z |
| 0.395 (0.125, 0.668) | |
ale_index_z |
| 0.630 (-0.340, 1.580) | |
canbics_index_z |
| -0.284 (-1.240, 0.683) | |
population_100_c |
| 0.141 (0.067, 0.217) | |
Hyper Parameters | |||
Precision for ui | 1.060 (0.563, 1.860) | 1.070 (0.551, 1.910) | 1.850 (0.722, 4.150) |
Phi for ui | 0.209 (0.034, 0.528) | 0.519 (0.160, 0.871) | 0.619 (0.170, 0.950) |
Pedestrians
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -1.460 (-1.890, -1.090) | -1.850 (-2.350, -1.390) | -1.580 (-2.170, -1.040) |
vandix_z | 0.654 (0.343, 0.970) | 0.693 (0.377, 1.010) | 0.665 (0.362, 0.966) |
ln_roads_km_c | 0.620 (0.269, 0.980) | 0.248 (-0.120, 0.622) | |
roads_prop_highway_arterial_z |
| 0.389 (0.148, 0.633) | |
ale_index_z |
| 1.570 (0.675, 2.450) | |
canbics_index_z |
| -0.983 (-1.880, -0.103) | |
population_100_c |
| 0.111 (0.039, 0.184) | |
Hyper Parameters | |||
Precision for ui | 0.626 (0.404, 0.925) | 0.676 (0.426, 1.020) | 1.080 (0.602, 1.820) |
Phi for ui | 0.267 (0.075, 0.549) | 0.474 (0.190, 0.778) | 0.373 (0.064, 0.783) |
Southeast
All
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | 1.270 (1.050, 1.480) | 0.651 (0.321, 0.973) | 1.560 (0.776, 2.340) |
vandix_z | 0.352 (0.057, 0.648) | 0.386 (0.111, 0.664) | 0.350 (0.094, 0.609) |
ln_roads_km_c | 0.779 (0.474, 1.090) | 0.490 (0.142, 0.844) | |
roads_prop_highway_arterial_z |
| 0.627 (0.331, 0.921) | |
ale_index_z |
| 1.630 (0.032, 3.250) | |
canbics_index_z |
| -0.712 (-2.390, 0.955) | |
population_100_c |
| 0.152 (0.040, 0.265) | |
Hyper Parameters | |||
Precision for ui | 0.750 (0.528, 1.030) | 0.827 (0.571, 1.150) | 0.994 (0.677, 1.400) |
Phi for ui | 0.396 (0.151, 0.689) | 0.691 (0.377, 0.923) | 0.771 (0.457, 0.961) |
Cyclists
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -2.000 (-2.780, -1.390) | -2.450 (-3.360, -1.640) | -0.882 (-2.090, 0.301) |
vandix_z | 0.597 (0.156, 1.060) | 0.755 (0.271, 1.250) | 0.828 (0.337, 1.330) |
ln_roads_km_c | 0.468 (-0.008, 0.948) | 0.202 (-0.535, 0.939) | |
roads_prop_highway_arterial_z |
| 0.331 (-0.323, 0.989) | |
ale_index_z |
| 1.940 (-0.715, 4.620) | |
canbics_index_z |
| 0.074 (-2.440, 2.590) | |
population_100_c |
| 0.288 (0.105, 0.473) | |
Hyper Parameters | |||
Precision for ui | 11.100 (0.759, 67.100) | 60.400 (0.709, 408.000) | 5,990.000 (0.933, 22,800.000) |
Phi for ui | 0.277 (0.018, 0.796) | 0.357 (0.028, 0.880) | 0.353 (0.018, 0.906) |
Pedestrians
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -1.870 (-2.560, -1.290) | -2.100 (-2.900, -1.400) | -1.380 (-3.000, 0.125) |
vandix_z | 0.545 (0.075, 1.020) | 0.597 (0.113, 1.080) | 0.606 (0.111, 1.110) |
ln_roads_km_c | 0.295 (-0.211, 0.801) | -0.070 (-0.840, 0.685) | |
roads_prop_highway_arterial_z |
| 0.307 (-0.349, 0.958) | |
ale_index_z |
| 0.615 (-2.210, 3.390) | |
canbics_index_z |
| -0.167 (-2.980, 2.700) | |
population_100_c |
| 0.233 (0.010, 0.464) | |
Hyper Parameters | |||
Precision for ui | 0.728 (0.368, 1.310) | 0.753 (0.372, 1.380) | 0.736 (0.349, 1.390) |
Phi for ui | 0.165 (0.011, 0.546) | 0.204 (0.014, 0.634) | 0.227 (0.016, 0.679) |
Northwest
All
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | 0.682 (0.281, 1.050) | 0.162 (-0.275, 0.573) | -5.270 (-21.900, 11.400) |
vandix_z | 0.281 (0.042, 0.521) | 0.362 (0.138, 0.591) | 0.222 (0.013, 0.433) |
ln_roads_km_c | 0.780 (0.469, 1.100) | 0.423 (-0.130, 0.986) | |
roads_prop_highway_arterial_z |
| 0.228 (-0.209, 0.663) | |
ale_index_z |
| 3.510 (1.690, 5.290) | |
canbics_index_z |
| -10.800 (-30.200, 8.550) | |
population_100_c |
| 0.372 (0.130, 0.614) | |
Hyper Parameters | |||
Precision for ui | 0.632 (0.381, 0.999) | 0.559 (0.323, 0.885) | 1.160 (0.688, 1.850) |
Phi for ui | 0.367 (0.054, 0.785) | 0.810 (0.430, 0.983) | 0.192 (0.002, 0.776) |
Cyclists
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -2.580 (-3.500, -1.670) | -2.630 (-3.650, -1.610) | -56.300 (-97.500, 6.250) |
vandix_z | 0.237 (-0.173, 0.649) | 0.246 (-0.180, 0.672) | 0.069 (-0.498, 1.250) |
ln_roads_km_c | -0.063 (-0.765, 0.638) | -0.993 (-5.480, 2.550) | |
roads_prop_highway_arterial_z |
| 0.654 (-1.950, 2.010) | |
ale_index_z |
| 2.980 (-7.640, 8.270) | |
canbics_index_z |
| -63.900 (-114.000, 6.580) | |
population_100_c |
| 0.758 (-0.010, 2.550) | |
Hyper Parameters | |||
Precision for ui | 1,270.000 (2.360, 7,980.000) | 1,280.000 (2.330, 8,040.000) | 1,380.000 (2.440, 8,570.000) |
Phi for ui | 0.343 (0.009, 0.927) | 0.343 (0.009, 0.927) | 0.344 (0.009, 0.927) |
Pedestrians
term | Unadjusted | Adjusted1 | Adjusted2 |
|---|---|---|---|
Fixed Effects | |||
(Intercept) | -1.470 (-2.210, -0.848) | -1.680 (-2.540, -0.951) | -5.510 (-23.500, 12.400) |
vandix_z | 0.275 (-0.006, 0.563) | 0.303 (0.010, 0.606) | 0.175 (-0.126, 0.478) |
ln_roads_km_c | 0.252 (-0.212, 0.740) | 0.143 (-0.699, 1.040) | |
roads_prop_highway_arterial_z |
| 0.258 (-0.452, 0.941) | |
ale_index_z |
| 4.090 (1.590, 6.640) | |
canbics_index_z |
| -9.180 (-30.100, 11.700) | |
population_100_c |
| 0.349 (0.032, 0.674) | |
Hyper Parameters | |||
Precision for ui | 1.540 (0.516, 3.900) | 1.370 (0.461, 3.420) | 2.330 (0.617, 6.960) |
Phi for ui | 0.341 (0.027, 0.849) | 0.450 (0.050, 0.917) | 0.305 (0.024, 0.804) |
end_time <- Sys.time()
end_time - start_time## Time difference of 23.88309 mins